DUDLEY KNO.X LIBRARY 
NAVAL POSTGRADUATE SCHOOL 
MONTEREY CA 93943-5101 



Approved for public release; distribution is unlimited. 

A Computaional and Experimental Investigation 
of the Propulsive and Lifting Characteristics of Oscillating Airfoils and Airfoil Combinations 

in Incompressible Flow 

by 

Kerry S. J'leace 

Lieutenant, United States Navy 
B.S., Georgia Institute of Technology, 1985 

Submitted in partial fulfillment 
of the requirements for the degree of 

AERONAUTICAL AND ASTRONAUTICAL ENGINEER 

from the 



NAVAL POSTGRADUATE SCHOOL 
September, 1992 



SECURITY CLASSIFICATION OF THIS PAGE 



REPORT DOCUMENTATION PAGE 


la. REPORT SECURITY CLASSIFICATION 
UNCLASSIFIED 


1b RESTRICTIVE MARKINGS 


2a. SECURITY CLASSIFICATION AUTHORITY 


3 DISTRIBUTION/AVAILABILITY OF REPORT 
Approved for public release; distribution is unlimited. 


2b. DECLASSIFICATION/DOWNGRADING SCHEDULE 


4. PERFORMING ORGANIZATION REPORT NUMBER(S) 


5 MONITORING ORGANIZATION REPORT NUMBER(S) 


6a NAME OF PERFORMING ORGANIZATION 
Naval Postgraduate School 


6b. OFFICE SYMBOL 
(If applicable) 
AA/PL 


7a NAME OF MONITORING ORGANIZATION 
Naval Postgraduate School 


6c. ADDRESS (City, State, and ZIP Code) 
Monterey, CA 93943-5000 


7b ADDRESS (City, State, and ZIP Code) 
Monterey, CA 93943-5000 


8a NAME OF FUNDING/SPONSORING 
ORGANIZATION 


8b. OFFICE SYMBOL 
(If applicable) 


9 PROCUREMENT INSTRUMENT IDENTIFICATION NUMBER 


8c ADDRESS {City, State, and ZIP Code) 


10. SOURCE OF FUNDING NUMBERS 


Program Element No 


Project No 


Task No 


Work Unit Accession 
Number 



1 1 . TITLE (Include Security Classification) 

A COMPUTATIONAL AND EXPERIMENTAL INVESTIGATION OF THE PROPULSIVE AND LIFTING CHARACTERISTICS OF 
OSCILLATING AIRFOILS AND AIRFOIL COMBINATIONS IN INCOMPRESSIBLE FLOW 



12. PERSONAL AUTHOR(S) Kerry S. Neace 



13a. TYPE OF REPORT 


13b. TIME COVERED 


14 DATE OF REPORT (year, month, day) 


15 PAGE COUNT 


Master’s Thesis 


From To 


1992, September 


123 



16 SUPPLEMENTARY NOTATION 

The views expressed in this thesis are those of the author and do not reflect the official policy or position of the Department of Defense or the U.S. 
Ilovernment. 



17 COSATI CODES 


18 SUBJECT TERMS (continue on reverse if necessary and identify by block number) 


FIELD 


GROUP 


SUBGROUP 


Phase relationships for oscillating airfoils. Propulsive forces (Katzmayr effect). 








Enhanced lift. Quasi-steady theory, Flow visualization 











19 ABSTRACT (continue on reverse if necessary and identify by block number) 

Computational and experimental methods have been used to systematically study one and two airfoils undergoing unsteady motion. First, a 
single airfoil analysis was done with the modified computer code, U2DIIF. Thrust, efficiency, and phase relationships were computed and 
:ompared to existing theoretical results. Furthermore, to help understand the dynamic stall process, relationships were developed between 
steady and quasi-steady pressure distributions for an airfoil undergoing a ramp motion. Next, an unsteady analysis for two airfoils was done with 
the modified computer code, USPOTF2. Again, thrust and efficienies for interfering, harmonically oscillating airfoils were computed and 
:ompared to existing theoretical results. Furthermore, an analysis was completed on the effects of a harmonically oscillating airfoil on the 
pressure gradient of a stationary airfoil. Finally, flow visualization experiments were conducted using a low speed smoke tunnel at the Naval 
Postgraduate School (NPS). This experiment demonstrated the effects of a thrust producing, oscillating airfoil on the formation of the wake 
cortices. Furthermore, a flow visualization experiment was conducted in the NPS low speed wind tunnel, which demonstrated the beneficial 
influence of a secondary airfoil oscillating in the vicinity of a stationary airfoil at high angle-of-attack. 



10. DISTRIBUTION/AVAILABILITY OF ABSTRACT 

QsUNCmSSIFIED/UNUMITED Q SAME AS REPORT 


] OTIC USERS 


21. ABSTRACT SECURITY CLASSIFICATION 
UNCLASSIFIED 


!2a. NAME OF RESPONSIBLE INDIVIDUAL 
Max F. Platzer 


22b. TELEPHONE (Include Area code) 
(408)646-2058 


22c. OFFICE SYMBOL 
AA/PL 



)D FORM 1473, 84 MAR 83 APR edition may be used until exhausted SECURITY CLASSIFICATION OF THIS PAGE 

All other editions are obsolete UNCLASSIFIED 



1 



ABSTRACT 



Computational and experimental methods have been used to 
systematically study one and two airfoils undergoing unsteady motion. First, 
a single airfoil analysis was done with the modified computer code, U2DIIF. 
Thrust, efficiency, and phase relationships were computed and compared to 
existing theoretical results. Furthermore, to help understand the dynamic stall 
process, relationships were developed between steady and quasi-steady 
pressure distributions for an airfoil undergoing a ramp motion. Next, an 
unsteady analysis for two airfoils was done with the modified computer code 
USPOTF2. Again, thrust and efficiencies for interfering, harmonically oscillating 
airfoils were computed and compared to existing theoretical results. 
Furthermore, an analysis was completed on the effects of a harmonically 
oscillating airfoil on the pressure gradient of a stationary airfoil. Finally, flow 
visualization experiments were conducted using a low speed smoke tunnel at 
the Navel Postgraduate School (NPS). This experiment demonstrated the 
effects of a thrust producing, oscillating airfoil on the formation of the wake 
vortices. Furthermore, a flow visualization experiment was conducted in the 
NPS low speed wind tunnel, which demonstrated the beneficial influence of a 
secondary airfoil oscillating in the vicinity of a stationary airfoil at high angle-of- 
attack. 



iii 



// 

t.l 

TABLE OF CONTENTS 

I. INTRODUCTION 1 

A. GENERAL 1 

B. SCOPE 2 

II. SINGLE AIRFOIL ANALYSIS 4 

A. HARMONIC MOTION 4 

1. Introduction 4 

2. Panel Code U2DIIF 7 

B. PROPULSIVE EFFICIENCY 16 

1. Introduction 16 

2. Theory 17 

3. Comparison to Flat Plate Theory 21 

4. Power Extraction for Two Degrees of Freedom 25 

C. UNSTEADY PRESSURE DISTRIBUTIONS 30 

1. Introduction 30 

2. Theory 30 

3. Description 32 

4. Results 33 



IV 



SCHOOL 



III. TWO AIRFOIL ANALYSIS 42 

A. COMPUTER CODE USPOTF2 42 

1. Modifications (USPOTF2A) 43 

a. Program Output 43 

b. Program Corrections 45 

c. Airfoil Motion 46 

d. Limitations 52 

B. COMPARISON TO EXISTING COMPUTATIONAL RESULTS . . 54 

C. PROPULSIVE EFFICIENCY 60 

1 . Introduction 60 

2. Theory 61 

3. Comparison to Flat Plate Theory 62 

D. UPWIND INFLUENCE STUDY 67 

1. Introduction 67 

2. Description 67 

3. Results 68 

4. Conclusion 69 



1. ai tY KNOX LIBRARY 
iVAVAL POSTGRADUATE 
MONTEREY CA 93943 -5 



IV. PROPULSIVE FLOW VISUALIZATION EXPERIMENT 73 

A. INTRODUCTION 73 

B. THEORY 74 

C. EXPERIMENTAL SETUP 76 



v 



1. Wave Propeller 76 

2. Wind Tunnel 77 

D. TEST PROCEDURE 79 

E. RESULTS AND DISCUSSION 80 

V. ENHANCED LIFT FLOW VISUALIZATION EXPERIMENT 89 

A. INTRODUCTION 89 

B. EXPERIMENTAL SETUP 90 

1. Wave Propeller and Stationary Wing 90 

2. Wind Tunnel 91 

3. Smoke Generator 91 

C. TEST PROCEDURE 92 

D. RESULTS AND DISCUSSION 92 

VI. CONCLUSION AND RECOMMENDATIONS 96 

A. SINGLE AIRFOIL ANALYSIS 96 

B. TWO AIRFOIL ANALYSIS 96 

C. FLOW VISUALIZATION EXPERIMENTS 97 

APPENDIX A 99 

APPENDIX B 104 

vi 



LIST OF REFERENCES 110 

INITIAL DISTRIBUTION LIST 112 



vii 



TABLE OF SYMBOLS 



b, c, I chord or 1/2 chord, as defined 
C L lift coefficient 

C Lo lift curve slope 

C D drag coefficient 

D time averaged drag over one cycle of motion 

D drag force 

h 0 plunge amplitude 

i denotes complex number: sqrt(-l) 

l L imaginary part of lift coefficient 

l M imaginary part of moment coefficient 

lm(x) imaginary part of x 

k reduced frequency 

L lift force 

L| imaginary part of lift 

L r real part of lift 

M moment 

N normal force 

P power 

q dynamic pressure 



VIII 



Re(x) 

Rl 

Rm 

T s 

» 

t 

u, v 
w 

a, AOA 

o 0 

A a 

n 

p 

0 

UJ 



real part of x 
real part of lift coefficient 
real part of moment coefficient 
suction force 
nondimensional time 
freestream velocity 

time averaged work over one cycle of motion 

angle-of-attack 

pitch amplitude 

change in angle-of-attack 

efficiency 

density 

phase angle between force and motion 
frequency of harmonic oscillation 



IX 



ACKNOWLEDGEMENTS 



The research for this thesis was conducted using the facilities of the 
Department of Aeronautics and Astronautics at the Naval Postgraduate School. 

I would like to give my sincere appreciation to professors M. F. Platzer and S. 
K. Hebbar, my thesis advisor and co-advisor, for their guidance, encouragement 
and many hours of council that led to the completion of this work. 

I would also like to thank the professional staff of the Department of 
Aeronautics and Astronautics, in particular Mr. Tony Cricelli for helping me 
through his gifted knowledge of the different computers in the department. 
Furthermore, the timely and precise technical assistance provided by Mr. John 
Molton, Mr. Ronald Ramaker, and Mr. Don Meeks made the experimental 
research a success. 

Finally, I would like to thank my lovely wife, Patricia, for her 
encouragement and support during my entire time here at the Naval 
Postgraduate School. 



x 



I. INTRODUCTION 



A. GENERAL 

In this thesis, several numerical methods were used to analyze the flow 
about one and two airfoils performing unsteady motion in an inviscid, 
incompressible fluid. First, the unsteady motion for a single airfoil was studied 
using the unsteady panel code, U2DIIF, by Teng [Ref. 1]. The primary purpose 
was to verify the phase lag relationships between the airfoil's motion and the 
build up of aerodynamic forces. To accomplish this, the time dependent output 
of the panel code was converted to harmonic output using a curve-fit algorithm. 
Furthermore, an extensive study on the production of thrust and associated 
efficiency for oscillating airfoils was completed using the U2DIIF code. Finally, 
a theory was presented that related quasi-steady ramp motion to a purely 
steady state phenomenon. 

Next, the vortex interaction between two airfoils was studied using a 
computer code, USPOTF2, developed by Pang [Ref. 12] for unsteady 
incompressible flow. Extensive modifications were made to increase the 
program capabilities. Comparisons were made to flat plate, linear theory using 
the modified code. Again, an analysis of the propulsive forces and efficiencies 
associated with two airfoils was completed using the modified code and a 
similar curve-fit algorithm developed to convert the time dependent output to 
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harmonic output. This led up to a systematic study involving the influence of 
an oscillating airfoil in the vicinity of a stationary airfoil. 

Finally, smoke flow visualization experiments were conducted using the 
available facilities located at the Naval Postgraduate School. These 
experiments were undertaken in an attempt to verify the two-dimensional 
theory presented here on the production of thrust associated with an oscillating 
airfoil. Furthermore, a visualization experiment was conducted to better 
understand the influence of an oscillating airfoil in the vicinity of a stationary 
airfoil at high angle-of-attack. 

B. SCOPE 

Chapter II contains the complete single airfoil, theoretical analysis using 
the panel code, U2DIIF. This chapter begins with the development of harmonic 
motion, includes an extensive verification with existing two-dimensional, linear 
theory, and concludes with a presented quasi-steady theory. Chapter III 
describes the complete two-airfoil, theoretical analysis. The development of 
the modified version of USPOTF2A is shown along with a complete description 
of the changes. Again, this code is verified with existing two-dimensional, 
linear theory, and used in a systematic study of interacting airfoils. Chapter IV 
and V incorporate some details of the flow visualization experiments and 
selected flow visuzalization photos of the results for the thrust and enhanced 
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lift investigation, respectively. Chapter VI contains the conclusion and 
recommendations. 
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II. SINGLE AIRFOIL ANALYSIS 



A. HARMONIC MOTION 
1 . Introduction 

The study of two dimensional, unsteady, harmonic motion involves 
an airfoil undergoing sinusoidal pitch or plunge oscillations. As the airfoil 
performs this oscillation, a complex flow field develops as shown in Figure 2.1 . 
As the airfoil increases its angle-of-attack, the pressure field around the airfoil 
changes. These changes create a disturbance in the boundary layer due to the 
presence of viscosity. This disturbance builds up on the airfoil surface as it 
travels to the trailing edge in the form of a vortex. The vortex then sheds into 
the medium with a circulation strength equal in magnitude to the increase in 
circulation about the airfoil, but opposite in direction. These disturbances are 
stored in the fluid because the shed vortices convect downstream at the local 
flow field velocity. The counter-rotating vortices induce a sinusoidal flow field 
which further changes the net lift. The result of this complex flow field is a 
time difference or delay in the airfoil's motion and the induced aerodynamic 
forces. This time difference is known as the phase lag, ip [Ref. 1]. 
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Figure 2.1 Oscillating Flow Field 

To simplify the calculations for this type of motion it is common to describe the 
airfoil position and the associated aerodynamic forces with complex variables. 
For pure plunge oscillations the vertical motion of the airfoil is described by the 
real part of the following equation: 



where, h 0 is a complex number, and a> is the frequency of oscillation. 
Similarly, the lift or moment is described by: 



where, L 0 is the theoretical quasi-steady lift given by the expression: 



h(t) = h 0 e^ 



( 2 . 1 ) 



L = L^re* 



(2.2) 




(2.3) 
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This is termed quasi-steady because the angle of attack is represented by h/U. 
The values of r and ip represent the magnitude and phase, respectively, of the 
true instantaneous lift relative to the quasi-steady lift. The variables r and tp in 
general depend on the reduced frequency k, the Mach number M, and the 
Reynolds number. For a nonviscous, incompressible fluid the values of r and 
tp will be functions only of k. A complete solution for the oscillating flat-plate 
airfoil in incompressible flow has been obtained by Kussner and Theodorsen and 
is reproduced from reference 3 in Figure 2.2. 
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Harmonic motion is usually characterized by non-dimensional 
variables. The amplitude for plunge motion is made non-dimensional by dividing 
by a reference length, usually the chord or half chord. Pitch amplitude is 
expressed in radians. The nondimensional frequency is termed the "reduced 
frequency" and is given by the following expression: 

k = (2.4) 

where, oj = frequency of oscillation, b= chord (or half-chord), U= free- 
stream velocity. 

2. Panel Code U2DIIF 

The computer code U2DIIF was developed by Teng [Ref. 2] for the 
study of unsteady, inviscid and incompressible flow over a single airfoil. The 
code was based on a technique called Panel Methods developed by Hess and 
Smith [Ref. 4] for steady potential flow problems. This method involves dividing 
the airfoil into many segments or panels. A uniform source and vorticity 
distribution is placed on each panel. In this code, the source strength is 
allowed to vary from panel to panel while the vorticity strength is held 
constant. This method is based on the fact that the singularity distributions 
automatically satisfy Laplace's equation, which is the governing equation for 
inviscid, incompressible flow. Furthermore, since the superposition principle 
applies to the Laplace equation, one can build complicated flow fields by a 
combination of simple flows. The unsteady potential flow model is based on 
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a method by Basu and Hancock [Ref. 5]. This unsteady model is governed by 
the Helmholtz vortex theorem which requires that any change in circulation 
around the airfoil must be matched by an equal and opposite change in vorticity 
in the wake. This is known as the vortex shedding process, which is modeled 
in this code by a vortex panel that is shed into the wake at each time step. 
The introduction of a wake model creates a non-linear, unsteady flow problem 
which is solved through an iterative process. No attempt is made here to 
reproduce the work by Teng or to explain the operation of the U2DIIF code, but 
instead the reader is encouraged to review reference 2. One limitation in the 
U2DIIF code which should be noted is the sensitivity to panel density and airfoil 
thickness. It can be shown that for thin airfoils of less than 8% thickness an 
increasing number of panels ( < 50) is required to accurately capture the leading 
edge suction peak. Furthermore, as the airfoil approaches a flat plate there is 
no amount of panels that will capture this peak. Although flat plates or 
NACA0001 airfoils are used in the following studies, it should be noted that the 
results of this code for airfoils of less than 5% thickness should be considered 
suspect unless it can be validated by other means. 

The code was verified by comparing the time dependent output for 
harmonic motion to experimental results obtained by Giesing [Ref 6.]. In order 
to further verify the U2DIIF code and to present the results for harmonic motion 
in a more useful manner, a program was written to convert the time dependent 
output of lift and moment histories to harmonic output using an iterative, curve- 
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fit algorithm (appendix A). Figure 2.3 shows an example of the time dependent 
output, the curve fit solution, and the associated phase angles between the 
motion and the aerodynamic forces. In order to convert the time dependent 
data to harmonic output each run must proceed for at least two cycles. The 
first cycle is discarded because of the transient effects associated with going 
from steady state to unsteady motion as seen in Figure 2.3. The second cycle 
of data (lift and moment histories) is then curve-fitted to the following 
expression: 

F{t) = Amp*Sin(ut + 4>) (2.5) 

where, Amp = amplitude of motion, cu = frequency of motion, and 0 = 
phase angle between the motion and the aerodynamic forces. 
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Figure 2.3 Sample Output. NACA0012 in pure pitch oscillations at a 
reduced frequency of 2, amplitude of 1 °, for 2 cycles. 
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3. Comparison to Flat Plate Theory 

Theoretical results for harmonic motion of a flat plate along with 
some experimental data were published by Halfman and reproduced in 
reference 7. The theoretical results for a flat plate oscillating in two degrees 
of freedom, h and a, were used to compare the U2DIIF code with harmonic 
output. The expressions of aerodynamic force and moment corresponding to 
the harmonic motion given in reference 7 are, respectively: 



Furthermore, the coordinate system used in the Halfman experiments was 
defined as positive lift acting down. A direct comparison of the U2DIIF 
harmonic output to these results was made by multiplying the U2DIIF output 
by 2 and adding 180° to the phase angle. This study was conducted with a 
panel density of 100 to accurately capture the suction peak for this very thin 
airfoil. The results are shown in Figures 2.4 through 2.7. It can be seen that 
there is very good agreement between the theoretical results and the U2DIIF 
code. 




( 2 . 6 ) 




(2.7) 
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Moment Coefficient VS. Reduced Frequency 
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Figure 2.5. Moment magnitude and phase comparison for a flat plate pitching 
about the .37c at an amplitude of 6.7°. 
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Figure 2.6. Lift magnitude and phase comparison for a flat plate in pure 
translation at a magnitude of .0833c. 
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translation at an amplitude of .0833c. 
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B. PROPULSIVE EFFICIENCY 



1 . Introduction 

Propulsion or thrust force from harmonically oscillating airfoils is a 
well known phenomenon. It has long been observed in nature from the low 
speed flapping of a bird's wing to the high speed flapping created by flying 
insects. The purpose of this study was to obtain the propulsive efficiencies 
associated with pure pitch and plunge harmonic motion using the nonlinear, 
unsteady panel code, U2DIIF. A comparison was made to the results obtained 
by Bosch [Ref. 8] for a flat plate undergoing harmonic motion. Bosch used a 
linear, analytical method to obtain the aerodynamic forces and propulsive 
efficiencies for a flat plate undergoing pure pitch and pure plunge harmonic 
motion. 

The propulsive forces associated with an oscillating foil have been 
experimentally measured in reference 18. In this reference, Scherer measured 
the propulsive forces and moments associated with a rigid foil of finite span 
undergoing large amplitude pitching, and plunging oscillations in water. This 
work was undertaken for the preliminary design of an oscillating foil propulsor 
suitable for use on a shallow-draft boat, such as a 'ski-barge.' The results are 
expressed in coefficient form, and may be used in a comparative analysis with 
the theory presented here using the appropriate three-dimensional corrections. 
Sufficient time did not permit a detailed comparison in this paper. 
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2. Theory 

In an inviscid, incompressible flow field propulsive forces develop over 
a harmonically oscillating airfoil. To better understand this process, Figure 2.8 
shows the forces that develop on a harmonically plunging airfoil. It can be seen 
that when the airfoil is moving vertically up or down the motion creates an 
induced velocity component. The force that develops remains perpendicular to 
the relative velocity. This allows for a forward or thrust component of force 
while the vertical component changes sign. This is an over-simplified 
explanation of a complex flow field. The actual flow field that develops will 
depend on the type of motion, the frequency of oscillation, and the intricate 
system of shedding vortices that store the kinetic energy of the motion in the 
fluid. 
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To understand how the mechanical energy of the motion is 
transferred to a propulsive force, it is necessary to develop the relationship for 
the required input work to the airfoil for both pitching and plunging motion. 
First, for a plunging airfoil the average work is defined as: 

r 

W = -( Re[L]*Re[h] dt (2.8) 

T o 



h = h 0 e , h = iu>h 0 e ht , T = — (2.9) 

Co 



The equation becomes, 




( 2 . 10 ) 

Re[(L R + iL,)][ coso) t + /sin o> f] Re[ iu h 0 (cosco t+is\nut)]dt 



After evaluation the equation reduces to, 



I Cl) d. 



( 2 . 11 ) 
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where, the imaginary lift and reduced frequency are given by, 



L, = pU 2 b lm(Ci) , 




( 2 . 12 ) 



Finally, the expression for average work reduces to, 



W = U 3 bklm(C J 



(2.13) 



This equation for average work agrees with reference 8 if we set the value 
h 0 /b = 1 .0. This equation shows the important functional relationships between 
the input variables of plunge amplitude, reduced frequency, and generated lift 
to the output variable of average work. Furthermore, it shows that the average 
work is proportional to the reduced frequency and the imaginary or 'out-of- 
phase' component of lift. 

The development for the relationship of average work for a pitching 
airfoil proceeds in a similar manner starting with the relation, 



T 




0 



(2.14) 



and the final result is, 
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W = a 0 p V*bklm{C M ) (2.15) 

This result again matches the equation in reference 8 if we set the value of 
<7 0 = 1.0 radian. This equation shows the functional relationship between the 
input variables of pitch amplitude, reduced frequency, and the imaginary 
component of moment with the output variable of average work. 

A measure of the efficiency with which work input is transferred into 
propulsive power is given by the following expression, 




W 



where, q = efficiency, D= average drag or propulsive force, v = velocity. 
This equation for efficiency is simply the drag (propulsive force) times the 
velocity divided by the expression for average work per unit time which is the 
average power output over average power required for the airfoil motion. It is 
clear the efficiency is primarily a function of the thrust that develops for a 
particular type of motion. Since the magnitude of thrust is proportional to the 
induced velocity ahead of the airfoil, the motion which provides the greatest 
induced velocity for the same work input will also provide the greatest thrust 
at the highest efficiency. The two types of motion studied in this paper are 
pitch, and plunge. Pitch motion induces a flow field in front of the airfoil that 
is proportional to the distance the pitch axis is from the airfoil. Furthermore, 
plunge motion can be thought of as pitch motion with the pitch axis at <». 
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This can be visualized as the motion of a fan blade with pitch motion similar to 
the blade hub, and plunge motion similar to the blade tip. Therefore, it should 
be no surprise that pure plunge motion provides more efficient propulsion than 
pure pitch. 

3. Comparison to Flat Plate Theory. 

Analytical results for a flat plate airfoil undergoing pure pitch and pure 
plunge motion for both propulsive force and efficiency were computed by 
Bosch in reference 8. To compare these results with the nonlinear code, 
U2DIIF, it was necessary to convert the efficiencies to aerodynamic forces. 
Furthermore, Bosch uses a reference length of b/2 (1/2 chord length) where the 
U2DIIF code uses b (chord length) as the reference length. Also, Bosch 
defines the average drag as follows, 

D=plv 2 C~ D (2.17) 

where 1 = 1/2 chord. Combining equations 13, 15, and 17 into 16 gives the 
foilwing relations: 

Cp _ 2C c ^2 -|gj 

" klm(CJ ’ Ap,un " " kim(Cj 

These equations for efficiency were used to compare the results obtained by 
Bosch. The U2DIIF code was executed at twice the reduced frequency as the 
one desired for comparison due to the different reference lengths. 
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In this study there were several problems associated with comparing 
analytical, linear results using a nonlinear code. First, the amplitudes used by 
Bosch: h/b = 1 .0, and a 0 = 1 .0 radian, (used for convenience) were well out of 
the linear range. For this reason much smaller amplitudes were used for the 
U2DIIF code: h/b = .05, and a 0 = .0873 radian, and the results were scaled 
appropriately. Scaling lift and moment results followed a linear relationship 
which implies that scaling the amplitude of motion by a factor results in scaling 
the forces up by the same factor. Drag follows a nonlinear relationship since 
drag is proportional to a 2 . Therefore, scaling the drag required multiplying the 
results from U2DIIF by the square of the scaling factor. After making all these 
conversions the U2DIIF harmonic output was compared to Bosch results in 
Figures 2.9 and 2.10 . 

Excellent agreement is shown in plunge for both the drag coefficient, 
and the propulsive efficiency. This is not the case for pitch. Good agreement 
is shown for propulsive efficiency in pitch, but the drag coefficient does not 
agree well at the higher reduced frequencies. The reason for this is most likely 
the result of scaling errors. When small amplitudes are used for pitch motion 
the resulting drag coefficient magnitudes are very small -.005. This is near 
the same magnitude as the accuracy of the code. Therefore, the inherent error 
in the code is scaled with the coefficient of drag magnitude. Plunge drag 
coefficient data are an order of magnitude larger than pitch and do not have the 
same scaling problem. 
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Figure 2.9 Average drag and efficiency for a NACA0009 airfoil undergoing 
pure plunge oscillations at an amplitude of .05c. The results were scaled up 
to 1.0c (xlOO) for comparison to Bosch. 
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Figure 2.10 Drag and efficiency for a NACA0009 airfoil undergoing pure 
pitch oscillations at an amplitude of 5.7°. The results were scaled up to an 
amplitude of 57 ° (xlOO) for comparison to Bosch. 



24 



4. Power Extraction for Two Degrees of Freedom 

When an airfoil undergoes both pitch and plunge motion it is 
oscillating with two degrees of freedom. For this condition it is possible for a 
phase relationship to exist between the two types of motion which results in 
a positive energy input into the airfoil. When this phase relationship exists 
energy is extracted from the airstream and transferred to the airfoil. This is the 
basic mechanism behind aerodynamic flutter. This may also be used as an 
efficient source of mechanical energy as described by McKinney and Delaurier; 
The Wingmill: An OscU/ating-Wing Windmill [Ref. 9]. The objective of this 
study was to compare the results obtained from the U2DIIF code for power 
extracted from airfoil motion in two degrees of freedom to the experimental 
results obtained in reference 9. First, the U2DIIF code was modified to perform 
the following harmonic motion, 

h = /7 0 sin (co /) , a = a 0 Sin(cof+4>) (2.19) 

where, 0 is the phase angle between plunging and pitching motions. Next, the 
equation for power extracted from the airstream is given by, 

P = f?(A/cos<x + (T s -D) sin a) + a M (2.20) 

where, N = normal force, D = drag, T 8 = suction force, and M = moment 
about the pitch axis. Assuming small perturbations and noting the opposite 
signs of T 8 and D, we can write, 
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( T s - D ) sina < A/cosa , 



Lift * N 



( 2 . 21 ) 



Therefore, the equation for power is reduced to, 

P = Lh * Ma (2.22) 

Equation 2.22 was used to obtain the power extracted during one 
cycle of motion. Experimental results are given in reference 9 for a wing 
undergoing the harmonic motion described in equation 2.19. The wing has a 
NACA0012 airfoil section and a rectangular planform of 20cm by 105cm. In 
the experiment h 0 was kept constant at 6 cm, and a 0 was set to 25°. The 
reduced frequency was given to be .361 based on the half chord and the phase 
angle was set to 90°. 

Figure 2.11 shows the airfoil motion for both angle-of-attack and 
translation as a function of time for one cycle. This plot shows the 90° phase 
difference between a and h. Figure 2.12 shows the aerodynamic forces in 
newtons and the power extracted in watts obtained from reference 9 and the 
U2DIIF code. Numerically, the experimental results differ from what we 
obtained by a factor of 2 both in the aerodynamic forces and the power 
extraction. They obtained a peak normal force coefficient equal to .897, while 
we obtained a peak lift coefficient of 2. This difference is due to the three- 
dimensional effects associated with the experiment combined with our 
approximation of the normal force coefficient. In reference 9 a lift-curve slope 
of C Lo = 4.30/rad was used to correct the two-dimensional Theodorsen theory. 
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If we apply the same correction to the theoretical lift-curve slope of c, 0 = 2/7 
used in U2DIIF along with the normal force correction, we obtain a normal force 
coefficient equal to 1.23. This is slightly greater then the results of reference 
9. This is what we expected considering the low Reynolds number and added 
mass term used in the wind tunnel experiments which would tend to lower the 
results from theoretical predictions. 

The power extraction curve followed the same trend as the 
aerodynamic forces, and again was off by a factor of 2. After applying the 
same corrections as above we obtain a peak power equal to 16.6 watts 
compared to the experimental value of 13.5 watts. It was encouraging to see 
the predictions from U2DIIF produce the same shape and phasing as the 
experimental results. Furthermore, we accurately predict positive work (power) 
over most of the cycle, as is expected for a phase angle of 90°. 
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Figure 2.11 Airfoil motion over two cycles for combined pitch and translation showing the 90 degree 
phasing between the two. 




10 




08 


3 


06 


i 


04 


S 

2 


02 


a 


0 


ti 




KJ 


-02 


*- 


< 


04 


£ 


•06 


'S 


-00 




•1 0 





POWER & AERODYNAMIC FORCES 




Nondimensional Time 



Figure 2.12 The aerodynamic forces (newtons) and power extracted (watts) 
for both the results obtained by McKinney [Ref. 9] and the results from the 
U2DIIF code respectively. 
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C. UNSTEADY PRESSURE DISTRIBUTIONS 



1 . Introduction 

The purpose of this study was to establish relationships between 
steady state pressure distributions and unsteady pressure distributions for a 
NACA001 2 airfoil using the panel code U2DIIF. The study centered on finding 
an unsteady pressure distribution on the upper surface of an airfoil undergoing 
a ramp motion which would match a pressure distribution from the same airfoil 
at some steady state value. Once the relationship was developed between 
steady and unsteady pressure distributions, the results could be used as an 
input in a direct boundary layer code. This would give a quick and efficient 
method of determining unsteady boundary layer profiles and the onset of 
dynamic stall using a simple panel code coupled with a direct boundary layer 
code. 

2. Theory 

The principal theory or hypothesis is that an airfoil undergoing a quasi- 
steady ramp motion can be analyzed as a steady-state process. The current 
thought process is that a ramp motion is a purely unsteady phenomenon which 
requires an unsteady analysis (ie. Navier Stokes) to determine the boundary 
layer and the onset of dynamic stall. In fact, as pointed out in reference 1 1, 
dynamic stall is strongly dependent upon airfoil geometry, Mach number, pitch 
rate, Reynolds number, state of the airfoil boundary layer, type of motion, etc. 
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Furthermore, the computation of unsteady, compressible flows can be a long 
and complicated process, even assuming the flow can be modeled 
appropriately. The next thought process is to assume that the boundary layer 
behaves in a quasi-steady fashion. This implies that the boundary layer reacts 
nearly instantaneously to the pressure distribution. More precisely, the 
boundary layer will react to the pressure gradients. The stronger the adverse 
pressure gradient the more likely the boundary layer tendency is to separate. 
To check this theory, a study was done in reference 10 where the unsteady 
pressure distribution for a particular motion was input into a steady direct 
boundary layer code. The direct boundary layer code can then accurately 
predict the onset of separation or flow reversal, although the code breaks down 
after separation due to the formulation of the direct boundary layer problem 
[Ref. 10]. The onset of separation for different pitch rates was compared to 
experimental results obtained by Chandrasekhara, Carr, and Ahmed at NASA 
Ames Research Center [Ref. 1 1]. This process covers all the variables stated 
above that influence dynamic stall with some limitations. First, the airfoil 
geometry, type of motion, and pitch rate are covered in the unsteady panel 
code. The effects of Reynolds number, and Mach number are covered in the 
direct boundary layer code because this code was modified to include the 
Prandtl-Glauert compressibility correction factor. The primary limitations in this 
analyis are that the three-dimensional effects are not considered, and at higher 
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pitch rates the flow becomes transonic making the compressibility correction 
less accurate. 

Once it has been established that the boundary layer behaves in a 
quasi-steady fashion for certain flows, a relationship can be obtained that links 
the quasi-steady pressure distribution on the upper surface to a purely steady 
state pressure distribution. For this to be true we must assume that the entire 
flow field is quasi-steady. This requires keeping the ramp motion sufficiently 
slow to satisfy this assumption. When an airfoil moves in a quasi-steady ramp 
motion the pressure field will require a finite amount of time to react as 
explained in the first section. For this reason the pressure distribution on the 
upper surface of the ramping airfoil will match a steady state pressure 
distribution at some earlier AOA. The purpose of this study was to verify this 
hypothesis and to determine the relationship between steady and quasi-steady 
flow fields. 

3. Description 

This study involved using the U2DIIF code for NACA0012 airfoils 
undergoing a ramp motion at various rise times to determine the relationship 
between steady and quasi-steady flows. It was first necessary to determine 
the range of rise times for which the assumption of quasi-steady flow would 
hold true. The rise time is defined as the time required for the airfoil to ramp 
from the initial AOA to the final AOA. The rise time is expressed in 
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nondimensional units by Teng [Ref 2.] and is defined as the time required for 
a fluid particle to travel from the leading edge to the trailing edge. 




c 



Since experimental results can normally be assumed quasi-steady due to 
mechanical limitations, the experimental results in reference 1 1 were used as 
a guide. In these experiments the ramp motion goes from 0 to 5° in 
0.001 67sec. The fluid particle takes 0.0005 sec to go from the leading edge 
to the trailing edge. Therefore, the nondimensional rise time is equal to 3.34. 
This is equivalent to a nondimensional pitch rate as defined in reference 1 1 of 
.027 units. Since this study involved ramp motions from 0 to 15°, rise times 
on the order of 10 units were considered quasi-steady. 

The study involved running the U2DIIF code for a particular rise time. 
The pressure distribution output was then compared to the steady state 
pressure distribution for AOA's of 6, 8, and 10 degrees. The comparison 
continued until the best match was found between the unsteady pressure 
distribution and the steady state distribution. In order to compare the pressure 
gradients on the upper surface, the pressure distribution from the U2DI IF output 
was differentiated using a four point central difference scheme. 

4. Results 

Figures 2.13 and 2.14 show a typical comparison study of the 
pressure distributions and pressure gradients. The steady state distribution is 
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plotted for a particular AOA along with several unsteady distributions for 
different AOA's until a 'best' match is found. In this particular study the match 
for a 10° AOA steady distribution was found at the unsteady AOA of 13.34° 
for a relatively fast rise time of 3 units. Figure 2.15 shows the AOA and the 
aerodynamic forces as a function of time for this ramp motion. This important 
result verifies that there does exist an unsteady pressure distribution that will 
match a steady state pressure distribution on the upper surface. 

One of the first results discovered in this study was that when the 
rise times were decreased from 10 to 3 units the best match became 
increasingly off. This was to be expected as the ramp motion increased to 
where the quasi-steady assumption would be invalid. It was found that rise 
times faster than 3 units would not match steady state pressure distributions. 
To show the trend, Figures 2.16 and 2.1 7 show the best match of a 10° AOA 
steady pressure distribution and the corresponding unsteady pressure 
distribution for a rise time of 10 and 3 units respectively. As can be seen, the 
match becomes progressively worse with the faster rise times. 

The numerical results for this study are tabulated in Table 2.1. 
Certain trends can be observed from these results. First, the steady state lift 
coefficient at a particular AOA is always greater than the corresponding 
unsteady lift coefficient at the same AOA. This is consistent with Theodorsen 
theory. Next, the peak unsteady lift coefficient is always greater than its 
corresponding steady state value and increases with increasing pitch rate. This 
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is consistent with the increased maximum AOA and corresponding lift obtained 
from dynamic motion over steady state results as shown in reference 1 1 . 

The AOA for which a particular match was found also followed 
certain trends. First, the delta AOA at the smaller steady state angle-of-attack 
(6 ° ) was almost independent of the rise time. This is equivalent to stating that 
it is independent of the pitch rate. A nominal value of two degrees was found 
between the unsteady and steady distributions for the smaller steady state 
AOA's. Next, the higher the steady state AOA, the more dependency on pitch 
rate was found. At 10° AOA steady the delta AOA went from 2.13 degrees 
for a rise time of 10 units to a delta AOA of 3.2 degrees at a rise time of 3 
units. 

When comparing this data to the experimental results of reference 1 1 , 
it was not straight forward. The best we could do was compare the effects 
of pitch rate on the onset of separation. Reference 11 notes that when the 
pitch rate is doubled from .03 to .05 units the onset in separation is delayed by 
2° (from 14° to 16°). Although we do not compute separation, we can 
compare the delay in matching the steady state values at an AOA of 14° for 
the pitch rates used in reference 1 1 . The results of this study are tabulated at 
the end of Table 2.1. As can be seen, we achieve the same delay of 2 ° when 
we double the pitch rate for a rise time of 14 to 7 units with a ramp change of 
20 °. 
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TABLE 2.1 Results of unsteady to steady pressure distribution matching. 



Steady 

AOA 


Steady 

c L 


Unsteady 

AOA 


Delta 

AOA 


Unsteady 

c L 


Rise Time =10 units (0° - 15") 


6 


.723 


7.72 


1.72 


.835 


8 


.960 


10.14 


2.14 


1.067 


10 


1.198 


12.13 


2.13 


1.254 


Rise Time = 8 units (O’ - 15°) 


6 


.723 


7.78 


1.78 


.856 


8 


.960 


10.12 


2.12 


1.07 


1C 


1.198 


12.33 


2.33 


1.26 


Rise Time = 6 units (O’ - 15") 


6 


.723 


7.87 


1.87 


.902 


8 


.960 


10.25 


2.25 


1.104 


10 


1.198 


12.51 


2.51 


1.281 


Rise Time = 4 units (0" - 15°) 


6 


.723 


7.9* 


1.9 


.999 


8 


.960 


10.36* 


2.36 


1.18 


10 


1.198 


12.82 


2.82 


1.32 


Rise Time = 3 units (0° - 15") 


6 


.723 


8.1** 


2.1 


1.1* 


8 


.960 


10.43* 


2.13 


1.26 


10 


1.198 


1.72 


3.2 


1.35 


Rise Time = 14 units (0° - 20") 


14 


1.87 


16.1 


2.1 




Rise Time = 7 units (0" - 20") 


14 


1.67 


18.1 


4.1 


1.8 
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Figure 2.13 Pressure distribution match for a steady AO A of 10 degrees and a rise time of 3 units. 
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Figure 2.14 Pressure gradient match for a steady AOA of 10 degrees and a rise Time of 3 units. 



Airfoil Motion 




Aerodynamic Forces 




Figure 2.15 Airfoil motion and aerodynamic forces for a NACA0012 airfoil 
undergoing a ramp motion from 0 to 15° in a rise time of 3 units. 



39 



co 

4—1 
• rH 

a 

d 

o 



<D 

6 



<u 

oo 

c 



O 

Pl| 

4-1 

d 

• rH 

X) 
r j 
u. 

o 



cx 

U 




VO 



ri 



o 



CO 

O 



vo 

o 



-+ 

o 



n 

o 



o 

♦ 

o 



<N 

I 



t 



VO 

I 



oo 

I 



xJ 

t-l 

o 

rd 

U 



X) 

c 

o 

£ 



<b- 



40 



Figure 2.16 Pressure distribution match for a steady AO A of 8 degrees and a rise time of 10 units. 
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Figure 2.17 Pressure distribution match for a steady AO A of 8 degrees and a rise time of 3 units. 



III. TWO AIRFOIL ANALYSIS 



A. COMPUTER CODE USPOTF2 

The computer code written by Pang [Ref. 12] was the primary method 
used in the following studies of two airfoils undergoing harmonic motion. The 
computer code, USPOTF2, was written to solve the potential flow for two 
airfoils executing unsteady motions in an inviscid, incompressible flow medium. 
This code is an extension of the single airfoil code, U2DIIF [Ref. 2], which uses 
the technique known as Panel Methods for steady flow and extends it to 
unsteady flow by introducing a wake model. This creates a non-linear problem 
due to the continuous shedding of vortices into the trailing wake. Furthermore, 
the presence of the second airfoil introduces a set of non-linear coupled 
equations for the Kutta condition. The solution requires an iterative procedure 
to compute the two vorticity strengths. Although no attempt is made here to 
reproduce the work by Pang, a general list of the modifications required to 
enhance the original code to a two-airfoil code is shown here. This will provide 
the reader with the necessary information to understand the modifications 
made by this author: 

1 . The establishment of five frames of reference: one fixed inertia frame of 
reference (global), two moving local frames of reference and two frozen 
local frames of reference. 
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2. Reformulation of the two Kutta conditions which are coupled non-linear. 

3. The creation of a new subroutine (NEWPOS) which transforms all 
coordinates in either of the two respective local frames of reference to 
the global frame of reference. 

4. The introduction of a more accurate method to obtain the velocity 
potential by integrating the velocity over smaller panels on the airfoil. 

5. Extension of the influence coefficient to include the effects of the 
second airfoil with its own peculiar wake. This also requires an 
introduction of an additional influence coefficient, that on the wake 
element due to the wake element from the other airfoil. 

6. The program is restricted to the following types of motion: 

• In-phase and out-of-phase Step Input 

• In-phase and out-of-phase Modified Ramp Input 

• In-phase and out-of-phase Translational Harmonic Oscillation 

• In-phase and out-of-phase Rotational Harmonic Oscillation 

• Sharp Edge Gust Field Penetration 
1 . Modifications (USPOTF2A) 

a. Program Output 

The code was originally written on the IBM mainframe at the 
Naval Postgraduate School. It has been transferred to the Stardent machine of 
the Department of Aeronautics and Astronautics. This allowed for increased 
storage capacity and a decrease in run times. For this reason, the total number 
of panels used to describe both airfoils has been increased from 200 to 400. 
The program was also transferred to the Iris workstation in the department. 
Although the program would compile on the workstation, this author could not 
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get the program to run on the Iris. It was believed to be a memory storage 
problem, but a thorough investigation was not completed. This program 
produces an extremely large amount of output to the screen. It is usually 
convenient on a Unix based machine to redirect screen output to a file during 
program runs. This file can get very large for long program runs, and may 
cause a problem in the temporary storage capacity. For this reason, the logical 
variable 'output' was added to the input file, FOROOI .DAT. When this variable 
is set to false, most of the screen output will not be printed, thus reducing the 
size of the output file. Furthermore, the important output has been separated 
into different files for easier analysis, and is not affected by the variable output. 
The following list describes the input/output files and the data they contain: 

1. FOR001.DAT: This is the input file (see Figure 3.1 ). 

2. FOR002.DAT: This is for user supplied airfoil coordinates, if 

desired. 

3. FOR003.DAT: This file contains the global coordinates of the first 
airfoil at each time step. 

4. FOR004.DAT: This file contains the global coordinates of the 
second airfoil at each time step. 

5. FOR007.DAT: This file contains the lift, drag, and moment 
coefficients for both airfoils at each time step. 

6. FOR008.DAT: This file contains the pressure coefficients for the 
first airfoil at each time step. 

7. FOR009.DAT: This file contains the pressure coefficients for the 
second airfoil at each time step. 
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8. FOR010.DAT: This file contains the first airfoil's core vortex (wake) 
positions at each time step. 

9. FOR011.DAT: This file contains the second airfoil's core vortex 
(wake) positions at each time step. 

10. F0R012.DAT: This file contains the AOA at each time step. 

11. F0R013.DAT: This file contains the DHY (translational motion) at 
each time step. 

12. F0R014.DAT: This file contains the required input for the Phase 
program. 

13. FOR020.DAT: This file contains the computed average pressure 
coefficient for the first airfoil. 

14. F0R021.DAT: This file contains the computed average pressure 
coefficient for the second airfoil. 

b. Program Corrections 

There were several errors noted in the original program and 
corrected. First, the code cannot be executed with only one airfoil as implied 
by the input variable NAIRF, or number of airfoils. Two airfoils must be defined 
to run the program. If single airfoil results are desired, one must position the 
second airfoil at a sufficient distance (approx. 30 chord lengths) away from the 
first airfoil to ensure no interference effects. This procedure was used to 

compare the USP0TF2 code against the single airfoil code, U2DIIF. At first, 

0 

the results did not match. This was when an error in the subroutine Press was 
discovered and corrected. Now, if one selects NGIES equal to 0.0, which is the 
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unsteady Kutta condition of equal pressures at the trailing edge panels, the two 
codes produce the same results. 

Next, an error in the convergence criteria was detected. 
Originally, the code would not reset the variable TOL, which is used by the 
code to set up the tolerance criteria for convergence. If, for a particular time 
step, the code required a higher tolerance to converge, the code would use this 
higher tolerance for follow-on time steps. This was corrected by resetting the 
variable TOL to its original value at the end of each time step. 

The original code computes the moment coefficient about the 
leading edge. Although this is not an error, it is usually desirable to have the 
moment coefficient about the pitch axis. This change was made in the 
modified code. 

c. Airfoil Motion 

This was the largest modification made to this code. As stated 
earlier, the original code was restricted to in-phase and out-of-phase motion. 
The code was modified so the airfoils can move independently of each other. 
This was done by adding the following variables into the code: DALP2, TCON2, 
FREQ2, and PIVOT2. The entire logic for airfoil motion was rewritten to include 
the new variables which allowed the airfoils to move independently. 
Corrections for this modification were made as needed in the subroutines. 
Figure 3.1 shows a sample input file, FOR001.DAT, required for the modified 
version of USPOTF2. Figures 3.2 and 3.3 demonstrate some of the new 
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capabilities of the modified code. In Figure 3.2, the first airfoil is stationary at 
10° AOA, while the second airfoil is oscillating in pure translation (DHY = .1) 
at a reduced frequency of 2.0. The aerodynamic forces clearly show an 
influence of the second airfoil's oscillation (drag coefficient) on the first airfoil's 
lift coefficient. In Figure 3.3, both airfoils are oscillating in pure pitch (DALP 
= 10°), but at different reduced frequencies. The top airfoil is oscillating at 
a reduced frequency of 4.0, while the bottom airfoil is oscillating at a reduced 
frequency of 1 .0. The aerodynamic forces are shown for 4 cycles of the top 
airfoil which is equivelant to 1 cycle of the bottom airfoil. The interaction of 
both airfoils can be seen in the aerodynamic forces with the stronger influence 
from the fast oscillating top airfoil. 

Another modification was made in the original code to allow for 
different size airfoils to be analyzed at the same time. The original code 
nondimensionalizes both airfoil's chord length to 1.0. The code was modified 
to include a new variable SCALE. This variable scales the second airfoil by the 
amount specified as a percent of the first airfoil. Figure 3.4 demonstrates the 
capability of this added feature. 



47 



Aug 17 1993 16:37:30 Stdin 


Page 1 


1 

2 

3 


NUMBER OF LINES FOR TITLE 
1 

TWO NACA00 1 2 OCSILLATING AIRFOILS 




5 


I FLAG 


NLOWER NUPPER 




6 


00, 


30, 30 




7 


NAIRFO 


, XSHIFT, YSHIFT, SCALE 




8 


2, 


1.2, .0, 1.0 




9 


NACA AIRFOIL TYPE, 




10 




12, 




11 




12, 




12 


ALP 1 , ALF 2 , DALPl , DALP2, TCCNl , TCON2 , 




13 


0.0, 0 


.0, 15.0, 0.0, 1.0, 0.0, 




14 


FREQ1 , 


FREQ2, PIVOT1, PIVOT2 
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OUTPUT 
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TRUE IF ONLY STEADY SOLUTION. FALSE OTHERWISE. 
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OUTPUT — 
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36 




X ( I ) , Y ( I ) COORDS. FOR BOTH AIRFOILS IN 




37 




FILE CODE 2. 
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40 




MAX AMPLITUDE OF AOA IN DEGREE FOR ROT. HARMON;'' MOTION. 
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FREQ1/2 
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46 


UGUST 
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48 


DELHXl /2 
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49 


DELHY1/2 
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(O TO BEGIN MOTION AT THE SAME TIME) . 
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TOL 


: TOLERANCE CRIERION FOR CONVERGENCE FOR (Uw) V and (Vw)k. 
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TAD J 


: FACTOR BY WHICH DTS WILL BE ADJUSTED. 




61 


SCLA 


: STEADY LIFT COEFF. FOR THE SINGLE AIRFOIL AT THF S r EC . AOA. 
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SCM 


: STEADY MOMENT COEFF. FOR THE SINGLE AIRFOIL. 
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SGAM 


: STEADY VORTICITY STRENGTH FOR THE SINGLE AIRFOIL. 
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NGIES 


: OPTION TO CHANGF THF UNSTEADY KUTTA CONDITION. 
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0 EQUAL PRESSURE AT THE TRAILING EDGE PANELS. 
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1 EQUAL TANGENTAL VELOCITIES AT THE TRAILING EDGE PANELS. 





Figure 3.1 Sample input file, FOROOI .DAT, required for the modified version 
of USP0TF2. 
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Wake Pattern 




Aerodynamic Forces 




Figure 3.2 Wake pattern and aerodynamic forces. AF1 -10 deg steady, AF2-0 
deg unsteady (DHY = . 1 ,k = 2) 
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Wake Pattern 
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Aerodynamic Forces 




Figure 3.3 Wake pattern and aerodynamic forces. AF1 -unsteady (DALP = 10 
deg, k = 1), AF2-unsteady (DALP = 10 deg, k = 4). 
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Figure 3.4 Wake pattern. Top airfoil scaled by .5 of bottom airfoil. Both undergoing a step inupt of 45 
degrees. 



d. Limitations 



There are several important limitations that must be understood 
thoroughly before attempting to use this program. First, the core vortices that 
make up the airfoil's wake must not intersect the panels that define the airfoil 
geometry. This situation is possible if one places the second airfoil in the wake 
of the first airfoil. The singularity nature of a core vortex intersecting a source 
can sometimes cause erroneous output. This can be noted in two ways. One 
way is to watch the screen output during the program run. If this situation 
arises, the code will not converge until the tolerance has been changed several 
times and is easily observed. Another way to identify this problem is to 
examine both the force, and moment data along with the wake patterns at the 
conclusion of the run. If the data or wake patterns are discontinuous, then 
there is a problem. This problem can be less severe if the panel density is 
increased and the time step is decreased. Unfortunately, if the wake is strong 
the interference effects are too great and the solution is suspect. Figure 3.5 
shows the wake pattern and aerodynamic forces of a stationary airfoil in the 
strong wake of an oscillating airfoil. This example clearly points out the biggest 
limitation with this program. 

Other limitations include a sensitivity to panel density, time step, 
and amplitude of motion. When adjusting these parameters one must watch 
the screen output to ensure the program can converge at each time step. 
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Main Title 




X - Axis 



Aerodynamic Forces 




Figure 3.5 Wake pattern and Aerodynamic forces showing vortex interference 
effects. AF1 unsteady, AF2 steady. 
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B. COMPARISON TO EXISTING COMPUTATIONAL RESULTS 



D.D. Liu and Z.X. Yao published results for vortex/wake flow studies in 
reference 1 3. For these results a numerical scheme was developed at Arizona 
State University that uses Panel Methods combined with an unsteady model 
that sheds vortices into the wake at each time step. The primary difference 
between this computer code and USPOTF2A is the treatment of the unsteady 
Kutta condition. Liu and Yao use a linearized model following the approach by 
Kim and Mook, where USPOTF2A uses the non-linear model developed by Basu 
and Hancock. 

Two case studies were chosen to compare the results. First, consider the 
case where two NACA0012 airfoils undergo a step change in AOA from 0 to 
5 ° . The second airfoil is placed .5c behind and .2c above the first. The results 
for wake patterns and aerodynamic forces for USPOTF2A are shown in Figure 
3.6. Liu and Yao have revised their results from reference 13. The current 
results are published here in Figure 3.7. The next case study involved the first 
airfoil undergoing a step change in AOA from 0 to 5° while the second airfoil 
remained stationary. Again, the second airfoil was placed .5c behind and .2c 
above the first. The results obtained by USPOTF2A are shown in Figure 3.8. 
The revised results obtained by Liu and Yao are shown in Figure 3.9. 

It is clear from the results that both codes follow the same trends. They 
do not produce the exact same results, as expected, due to the different 
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formulation of the problem. Both results essentially produce the same 
interference effect and approach steady state conditions. This example is a 
good verification for both numerical codes. 
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WAKE PATTERN 




LIFT COEFICIENT VS TIME 




Figure 3.6 USPOTF2A wake pattern and lift history output for two 
NACA0012 undergoing a step change in AOA from 0 to 5°. 
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Figure 3.7 Liu's output of lift histories for two NACAOOI 2 airfoils undergoing 
a step input from 0 to 5°. 
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Wake Pattern 




Aerodynamic Forces 




Figure 3.8 USPOTF2A wake pattern and lift history output for two 
NACA0012 airfoils. AF-1 undergoing a step change in AOA from 0 to 5°. 
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Figure 3.9 Liu's output of lift histories for two NACA001 2 airfoils, with the 
fist undergoing a step input from 0 to 5°. 



59 



C. PROPULSIVE EFFICIENCY 



1 . Introduction 

The purpose of this section is to develop relations for propulsive 
forces and efficiencies for two airfoils in tandem. There are two configurations 
discussed in this section (Fig. 3.10). One, the first airfoil oscillates while the 
second remains stationary. Two, the first airfoil remains stationary, while the 
second airfoil oscillates. The computational study will only examine the second 
configuration where the airfoil oscillates in either pitch or plunge. Due to the 
limitations of the code USPOTF2A, the case where the first airfoil oscillates is 
discussed with limited computations. The results of the numerical output were 
compared to the analytical results obtained by Bosch in reference 8. 
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Figure 3.10 The two configurations studied. 
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2. Theory 

The theory for the development of propulsive forces and efficiencies 
for two airfoils in tandem follows closely to the single airfoil case presented in 
Chapter II. The relationships for work done to the fluid by the two airfoils 
remain the same as for the single airfoil case since only one of the airfoils is in 
motion. The relationships for work done by plunge motion is given by equation 
2.15 and for pitch motion by equation 2.17. The primary difference between 
the single airfoil analysis and the two airfoil analysis is in the relationship for 
efficiency. Although only one airfoil is in motion, both airfoils are capable of 
producing thrust, and the expression for efficiency becomes: 

_ (5^5,) y (3.i) 

w 

where D 1( and D 2 are the drag or propulsive forces from the first and second 
airfoil, respectively. 

From the expression for efficiency in equation 3.1, the differences 
between the two configurations can be examined. In the first configuration, 
the pressure disturbances from the first airfoil are felt by the second airfoil 
downstream, or in the first airfoil's wake. This will produce a much stronger 
effect than the second configuration; where the pressure disturbances must 
travel upstream to have an effect on the stationary airfoil. Since the production 
of thrust by the stationary airfoil is entirely dependent on the influence of the 
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oscillating airfoil, it is clear that the first configuration will have a greater 
efficiency than the second configuration. Furthermore, as stated in chapter 2, 
plunge motion has a stronger influence than pitch motion; and the efficiencies 
associated with plunge motion are greater than those associated with pitch. 
Bosch concluded this in reference 8. Efficiencies for the first configuration near 
.9 were obtained for both types of motion with reduced frequencies between 
1 and 2. The second configuration reached an efficiency near .5 for plunge 
motion and .1 for pitch motion for the same range of reduced frequencies. 

3. Comparison to Flat Plate Theory 

Analytical results for two flat plates undergoing pure pitch and pure 
plunge motion for both configurations were computed by Bosch [Ref 8.]. To 
compare the results using USPOTF2A, the same conversions and scaling of the 
motion as for the single airfoil case were applied. Special consideration was 
given to the panel density and the time step to reduce the effects of the first 
airfoil's wake impinging upon the second airfoil. It was shown that when the 
first airfoil's wake is 'weak', it has a minimal effect on the solution or wake 
pattern of the second airfoil. The aerodynamic output was carefully observed 
proceeding each case study to ensure a smooth and continuous curve; this 
would confirm the impinging wake did not adversely effect the results. The 
time dependent output was then converted to harmonic output following similar 
procedures for the single airfoil case using the two airfoil phase-shift code 
presented in Appendex B. Then, the magnitudes of lift and moment were 
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divided into real and imaginary parts and plotted against reduced frequency. 
The results for pitch motion are shown in Figure 3.11; the first plot (C L ') 
represents the real part, and the second plot (C L ") represents the imaginary 
part. The agreement is particularly good at the low frequencies and begins to 
separate at the higher frequencies. 

The comparison for thrust and efficiency was done using the second 
configuration with plunge motion. The expressions for efficiency are given by 
the following relations: 

C D, * C D 1 _ 2(C 0| + C D; ) (3 2) 

r ' p " ch ' klm(CJ ' r ' pluns ‘ ‘ klm(CJ 

Figure 3.12 shows the propulsive force and efficiency comparisons. Here, the 
results match well over the entire frequency range with the computed 
efficiency always being less than that predicted by Bosch. This was 
anticipated because the non-linear effects in the USPOTF2A code would result 
in a lower efficiency. Furthermore, the consequence of plunge motion 
producing higher aerodynamic forces for the same frequency as pitch motion, 
results in less scaling errors. 

A final study was done to try and achieve the results by Bosch for the 
first configuration where the efficiencies are much higher. Unfortunately, as 
stated this code is unable to converge when the discrete vortices collide with 
the airfoil. Therefore, the second airfoil was placed 0.2c in front and 1c above 
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the stationary airfoil. In this configuration, the trailing vortices from the 
oscillating airfoil would not interfere with the stationary airfoil. The results 
were recorded for plunging motion (DHY = .05) at a reduced frequency of 4.0, 
and then the oscillating airfoil was brought incrementally closer to the 
stationary airfoil. This continued until the vortices interacted at 0.3c above the 
stationary airfoil. The efficiency went from /7 = .48, to n = - 6 at 1c and .3c 
above the stationary airfoil, respectively. This study shows the correct trend 
for the efficiency under the limitations of the code. 

Overall, the USPOTF2A program was shown to be an effective way 
to predict the aerodynamic forces, and efficiencies for two interacting airfoils. 
Unfortunately, the limitations of the program prevent some interesting case 
studies of strong wake interaction. Finally, one must be cautioned to observe 
the output with care when studying airfoils in close proximity to the wakes. 
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Figure 3.11 Lift coefficient as a function of the reduced frequency k. 
Pitching motion using configuration 2. 
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reduced frequency k. Plunging motion using configuration 2. 
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D. UPWIND INFLUENCE STUDY 



1 . Introduction 

It has been of interest to understand and study the effects of an 
oscillating airfoil in the vicinity of a stationary airfoil. The purpose of such a 
study is to see if there exists a configuration that would enhance the lift on the 
stationary airfoil. This is of primary interest in the post-stall region, where an 
oscillating airfoil may have the beneficial effect of delaying boundary layer 
separation and increasing the steady state lift. In theory, the oscillating airfoil 
will produce a certain amount of thrust, thereby increasing the flow velocity 
and creating a more favorable pressure gradient on the stationary airfoil. 
Although this study is done with the inviscid code, USPOTF2A, useful 
information on the pressure distributions and gradients are obtained from such 
an analysis. This information can be used to understand the degree of 
influence and trends in the integrated forces that determine the effects on a 
certain configuration. Due to the limitations of the code, only configurations 
where the oscillating airfoil was behind the stationary airfoil were considered. 

2. Description 

This study used the modified code, USPOTF2A, to determine the 
effects of an oscillating airfoil on a stationary airfoil at a high AOA. Two 
NACA0009 airfoils were used for this study with the first remaining stationary 
at 10° AOA. The second airfoil was set to oscillate in plunge (DHY = .2) at a 
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reduced frequency of 4.0. The amplitude of motion, DHY, was limited to .2 
because of the wake interference problem. The high reduced frequency was 
chosen to get the maximum effect with the relatively small amplitude of motion 
used in this study. The study involved changing both the position and mean 
AOA of the oscillating airfoil; then recording the change in pressure distribution 
and average lift of the stationary airfoil over steady state values. T wo positions 
for the oscillating airfoil were studied. First, the oscillating airfoil was placed 
slightly behind (XSHIFT = .1c) and above (YSHIFT = .6c) the stationary airfoil. 
Here, the mean AOA of the oscillating airfoil was changed from 0° to 10°. 
The other position considered was placing the oscillating airfoil slightly behind 
(XSHIFT = .1c) and below (YSHIFT = -.6c) the stationary airfoil. Here, the 
mean AOA of the oscillating airfoil varied from 0°, 5°, 10°, and 15°. 

3. Results 

The numerical results of this study are tabulated in table 3.1. In the 
first configuration, where the oscillating airfoil was located above the stationary 
airfoil, the influence was minimal. The increase in lift ranged from 7% to 10%. 
The important result was the influence of AOA. For a greater mean AOA on 
the oscillating airfoil, the delta lift on the oscillating airfoil was increased 
dramatically along with the influence on the stationary airfoil. Furthermore, the 
production of thrust decreased with increasing mean AOA. The results of the 
second configuration, where the oscillating airfoil was below the stationary 
airfoil, proved to be more effective. The increase in lift ranged from 10% to 
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20%. Again, the influence of mean AOA followed the same trends as the first 
configuration. 

A typical case study is displayed in Figures 3.13 and 3.14. Figure 
3.13 shows the steady state pressure distribution super-imposed on the 
average pressure distribution for one cycle of motion. Furthermore, the 
pressure gradients for both cases are displayed on the same graph. Figure 3.14 
shows the aerodynamic forces for two cycles of motion. The second cycle 
was analyzed for average data. It is interesting to note the cyclic behavior of 
the stationary airfoil's lift coefficient which is designated as AF-1 . Clearly, the 
influence of the oscillating airfoil on the stationary airfoil is shown in this graph. 

4. Conclusion 

Several important conclusions can be made from the above results. 
First, there is certainly a favorable influence from an oscillating airfoil in the 
vicinity of a stationary airfoil. Furthermore, the AOA and position of the 
oscillating airfoil have a strong influence on the effectiveness of the 
configuration. It was shown that the most effective configuration was having 
the oscillator at a high mean AOA and below the stationary airfoil. Finally, it 
can be concluded that the upwind influence can only produce a maximum of 
17% change in lift for an unstalled mean AOA of 10°. Unfortunately, the 
limitations of the program prevented studying many interesting configurations. 
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TABLE 3.1: Numerical results of oscillating airfoil study. 
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Figure 3.13 Pressure and gradient distribution for two NACA0009 airfoils. AF-1 is stationary at 10" 
AOA, while AF-2 is oscillating in plunge: A0A = 10\ DHY = .2, k=4.0, XSHIFT = .1, YSHIFT = -.6. 
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Figure 3.14 Aerodynamic forces for two NACA0009 airfoils. AF-1 is stationary at 10° AO A, while AF-2 
is oscillating in plunge: A0A=10°, DHY = .2, k = 4.0, XSHIFT = .1, YSHIFT = -.6. 



IV. PROPULSIVE FLOW VISUALIZATION EXPERIMENT 



A. INTRODUCTION 

The purpose of this experiment was to document the production of thrust 
by the wave propeller originally built at the Naval Postgraduate School by Carl 
Dane [Ref. 15]. This was a preliminary experiment to better understand the 
vortex pattern produced by a wave propeller, and to examine the production of 
thrust using smoke flow visualization techniques. 

An explanation of what constitutes a propulsive vortical signature along 
with smoke flow visualization of the propulsive vortical patterns is given in 
reference 16. In this reference, the explanation is given by contrasting the 
vortical pattern produced by a cylinder (drag) with the vortical pattern produced 
by a pitching airfoil (thrust). The cylinder produced a vortical street where the 
top row of vortices rotated clockwise and the bottom row of vortices rotated 
counterclockwise. This pattern induces a velocity component in the upstream 
direction (Biot-Savart law). In contrast, the pitching airfoil produced a 
counterclockwise rotating vortex street on the top row and a clockwise rotating 
vortex street on the bottom row. This pattern induces a velocity component 
in the downstream direction. Reproduction of the flow visualization data from 
reference 16 is shown in Figure 4.1. 
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Figure 4.1 Results from reference 16. Top vortex street shed from a circular 
cylinder (drag). Bottom vortex street generated by an airfoil pitching about the 
quarter-chord (thrust). 



B. THEORY 

A comparison was done using the incompressible panel code, U2DIIF. The 
purpose of this study was to examine the vortical pattern produced by the 
panel code, and determine if the vortical signature matched experimental 
results. The input to the panel code was set up to best match the conditions 
of the experiment described in the next section. The panel code was run using 
a plunge amplitude, h 0 /c, equal to .364, a reduced frequency of 2.5 and a 0° 
mean AOA. The results of the vortical pattern are shown in Figure 4.2. Aside 
from the starting vortex, this is clearly a thrust producing vortical street. 
Furthermore, the vortical pattern is similar to that produced by the experiment 
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Figure 4.2 Vortical pattern produced by the panel code, U2DIIF. h 0 /c = .364, and k = 2.5. 



shown in Figures 4. 1 0 and 4. 1 1 . Finally, the theory predicts a drag coefficient 
(thrust) of -.376. 

C. EXPERIMENTAL SETUP 
1 . Wave Propeller 

The wave propeller, originally proposed by Wilhelm Schmidt [Ref. 1 4], 
used in this paper was a modification of the original construction made by Carl 
Dane [Ref. 15]. The original mechanism was designed to perform a circular 
motion while holding the wing at a constant angle-of-attack. The primary 
modification was the construction of a new wing. The new wing was made 
from a NACA0012 airfoil section and consisted of a 5.5" chord and a 22" 
span. The wing was built from a foam core covered with a thin plywood skin 
and finished with a layer of glass-epoxy composite for added fatigue strength. 
The wave propeller's drive mechanism was rebuilt using added bearings and 
tighter fittings to allow for less binding at high rotational speeds. The final 
configuration for the wave propeller is shown in Figure 4.3. The drive motor 
used was a reversible DC motor rated at 24 volts, and 5 amps. The power 
supply used was regulated voltage DC power supply rated at 0-35 volts, and 
2.5 amps. The motor and power supply combination could provide the wave 
propeller with a maximum rotational speed of 1 500 revolutions per minute 
(rpm) when placed horizontal. 
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2. Wind Tunnel 



The wind tunnel used in this experiment was a very low speed, low 
turbulence smoke tunnel at the Naval Postgraduate School. This smoke tunnel 
is of indraft type and was designed to be a scaled model of an existing smoke 
tunnel at the NASA Ames Research Center. It is made of plexiglass walls and 
has a contraction ratio of 2.8:1. The motor provides wind tunnel velocities 
between 0 and 10 feet per second (fps). The smoke was created using a 
Rosco smoke generator and piped into the tunnel as a single stream for the 
streak line flow visualization experiment. The smoke was directed into a rake 
for the flow field visualization experiment. Figure 4.4 is a photograph of the 
wind tunnel and smoke rake used in this experiment. 
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Figure 4.3 Wave Propeller. 




Figure 4.4 Low speed smoke tunnel. 
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D. TEST PROCEDURE 



Testing was done in the low speed smoke tunnel for two different 
conditions. First, with the wind tunnel off (0 velocity), and second with the 
wind tunnel operating at 1 0 fps. In the first configuration with the wind tunnel 
off, the wave propeller was placed in the wind tunnel and the tunnel was filled 
with smoke; then, the wave propeller was turned on to its maximum rotational 
speed of 1 500 rpm and a rotation diameter of 2 inches. The purpose was to 
see if the wave propeller would draw the smoke through the tunnel like a fan, 
thus showing the production of thrust by the wave propeller. 

In the second tunnel condition, the smoke tunnel was turned on to its 
maximum power bringing the tunnel speed to 10 fps (chord Reynolds number 
of 26,500). The wave propeller was oscillated at three different rotational 
speeds: 165 rpm, 620 rpm, and 1085 rpm. This equated to reduced 

frequencies of .792, 2.97, and 5.2 respectively. Furthermore, the angle-of- 
attack of the wave propeller was varied between 0, 5, 1 0, and 20 degrees with 
a rotation diameter of 2 inches. For these configurations, smoke was used to 
visualize the flow field in two ways. First, smoke was introduced in one tube 
to visualize the streak line. Next, a rake was used to visualize the entire 2-D 
flow field. 

Photos were taken using a Nikon 35mm camera and Kodak TMAX-400 
ASA black and white film. The film speed was set to 1/500 seconds with an 
aperture setting of 2.0 for the light conditions used at low rotational velocities. 
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At higher rotational velocities film speeds of 1 /1 000 and 1 / 2000 seconds were 
used with aperture settings of 2.0 and 1 .4 respectively. 

E. RESULTS AND DISCUSSION 

The results for the first tunnel condition flow visualization experiment is 
shown in Figure 4.5. Here, the wind tunnel is turned off and the wave 
propeller is oscillating at approximately 1 500 rpm. It was shown that the wave 
propeller accelerated the smoke through the tunnel and out the exit similar to 
a fan blade. This experiment was done for the cases where the wave propeller 
oscillated both clockwise and counterclockwise with the same result. This was 
clear evidence the wave propeller produced thrust; hence, showing that circular 
motion can simulate plunge motion which is known to produce thrust. 

The results of the second tunnel condition flow visualization experiment 
are shown in Figures 4.6 - 4.20. Figure 4.6 shows the wave propeller at 
steady state and 0° AOA. It can be seen that the airfoil produces a highly 
turbulent wake at the very low Reynolds number (26,500) generated in this 
tunnel. The boundary layer would still remain attached, which was an 
improvement over the earlier airfoil design. Figures 4.7 through 4. 1 4 show the 
streak line flow visualization experiment at various AOA's, and frequencies. 
These photos show several important features about the street vortices 
produced by different oscillating conditions. Most of these pictures reveal the 
propulsive vortical street pattern discovered in reference 1 6. Unfortunately, the 
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wave propeller was constructed to oscillate at too high an amplitude (diameter) 
for this tunnel, and many of the vortical patterns would hit the walls and 
disperse. In Figure 4.10, the vortical pattern shows that the bottom vortex is 
rotating clockwise, and the top vortex is rotating counterclockwise, which is 
a thrust producing vortical street. It was observed that the AOA and frequency 
greatly affect the vortical strength (size). Increasing the AOA and frequency 
led to an increase in wake vorticity. Furthermore, these increases led to a 
stronger influence on the streak line in front of the wave propeller. Increasing 
the AOA beyond 5° led to stall, and subsequently, a decrease in thrust. Figure 
4.13 shows the wave propeller in a stalled condition (20° AOA), and the 
rotating vortices are both weaker and positioned such that they would not 
induce significant thrust. In fact, Figure 4.14 shows a dynamic stall vortex 
building up on the leading edge. 

Figure 4.15-4.18 contain the rake flow visualization experiment. Figure 
4.15 and 4.16 display the propulsive vortical pattern produced by the wave 
propeller for the 2-D flow field. Figures 4. 1 7 and 4. 1 8 indicate the large region 
of separation on the airfoil's upper surface. Again, the dynamic stall vortex is 
seen to form at the leading edge in Figure 4. 1 7. These pictures best reveal the 
stalled condition of the wave propeller operating at low Reynolds number and 
high amplitude. Both of these conditions were limitations in the design and test 
equipment. 
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Figure 4.5 Wave propeller operating at 1500 rpm, 0° AOA with the wind 
tunnel off. 







Figure 4.8 RPM = 165, AOA = 5°, rotating counterclockwise. 
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Figure 4.9 RPM = 620, AOA = 5°, rotating clockwise 
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Figure 4.11 RPM = 620, AOA = 10°, rotating counterclockwise 




Figure 4.12 RPM = 1085, AOA = 10°, rotating counterclockwise 
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Figure 4.13 RPM = 620, AOA = 20°, rotating counterclockwise. 




Figure 4.14 RPM = 1085, AOA = 20°, rotating counterclockwise. 
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Figure 4.17 RPM = 1085, AOA = 0°, rotating counterclockwise 




Figure 4.18 RPM = 1085, AOA = 0°, rotating counterclockwise 



88 





V. ENHANCED LIFT FLOW VISUALIZATION EXPERIMENT 



A. INTRODUCTION 

Experimental research in high angle-of-attack flight has long been an area 
of great interest. The area of active research centers around understanding and 
controlling the boundary layer to enhance lift in the post-stall region. To this 
end, many steady-state boundary layer control devices have been installed in 
an attempt to delay separation. 

Boundary layer control has also been attempted through unsteady 
excitation mechanisms. The most promising method is the wave propeller first 
suggested by Wilhelm Schmidt in the German Journal of Flight Sciences in 
1965 [Ref. 14]. Schmidt's wave propeller consisted of a single wing that 
performed a plunging motion perpendicular to the freestream. The oscillating 
wing was mounted between two stationary main lifting wings. In this 
arrangement the angle-of-attack of the wave propeller was varied to achieve 
optimum performance. Results from Schmidt's work indicate that with the 
optimum configuration the steady-state stall angle-of-attack could be increased 
to beyond 25° and the corresponding lift coefficient increased by four times. 

The results on flow visualization in this report was a continuation of the 
on-going research at the Naval Postgraduate School to better understand the 
beneficial effects of a wave propeller operating in the vicinity of a main wing. 
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The initial research and construction of the wave propeller was undertaken by 
a former NPS student, Carl Dane. Carl Dane's thesis [Ref. 15] contains the 
results of pressure distributions on a main wing with the wave propeller 
configured to oscillate aft and slightly above the wing. Unfortunately, the wave 
propeller was designed poorly due to construction limitations and no beneficial 
effects were observed. The purpose of this research was to build a more 
efficient wave propeller; then, through flow visualization, understand the 
effects of the wave propeller on a stationary wing. 

B. EXPERIMENTAL SETUP 

1 . Wave Propeller and Stationary Wing 

The wave propeller was the same as described in Chapter IV (Figure 
4.3). In this experiment, the wave propeller was mounted vertically in the NPS 
low speed wind tunnel behind a stationary wing. The wing had a 1 2-inch chord 
and a 28-inch span, the vertical length of the test section [Ref. 1 5]. The wing 
was constructed of a NACA 66(215)-216 airfoil section and could be set 
between 0 and 20 degrees AOA. Unfortunately, the wave propeller operated 
less efficiently in the vertical position. This was most likely a result of the 
increased weight on the bearings. A new power supply was used to provide 
increased rotational speeds of the wave propeller. A DC reversible power 
supply rated at a constant 48 volts, and 1 5 amps was used. With this power 
supply, the wave propeller operated at a rotational speed of 600 RPM. 
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2. Wind Tunnel 



The experiment was conducted in the 32x45 inch low speed wind 
tunnel located at the Navel Postgraduate School [Ref 17]. This tunnel was 
designed at the Aerolab Development Company of Pasadena, California as a 
closed circuit, single-return wind tunnel that is 64 feet long and 25.5 feet wide. 
The tunnel is powered by a 100 horse power motor which drives a three blade 
variable pitch fan by a four-speed International truck transmission. This motor 
and transmission can efficiently drive tunnel velocities up to 180 mph. Fan 
induced swirl is removed by eight stator blades located directly downstream of 
the fan. The tunnel cross section gradually expands from the fan to the settling 
chamber while turning through three sets of 90° corner vanes. The settling 
chamber includes two fine wire mesh anti-turbulence screens designed to break 
up large turbulent fluctuations. The flow accelerates to the test section though 
a 10:1 contraction cone. 

The test section is made up of transparent glass sidewalls (access doors) 
and upper wall for illumination and viewing. The test section is designed to 
operate at atmospheric pressure, and has a breather slot located aft of the test 
section to compensate for leakage losses throughout the tunnel. 

3. Smoke Generator 

The smoke was generated using the same Rosco smoke machine 
described in Chapter IV. The smoke was piped from the top of the tunnel into 
a blast tube. The blast tube directed the smoke into the tunnel slightly 
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upstream of the test section. In this configuration, the smoke was injected into 
the test section directly at the stationary wing. 

C. TEST PROCEDURE 

The wave propeller was mounted vertically from the top of the wind 
tunnel to extend down into the test section slightly behind the stationary wing 
(Figures 5.1 and 5.2). In this position, the wave propeller could rotate both 
clockwise or counterclockwise at various angle-of-attack settings. 

The tunnel was operated at speeds between 10-30 fps. The smoke was 
injected in the tunnel slightly upstream of the stationary wing. At this point, 
the wing angle-of-attack was increased until the boundary layer was just 
beginning to separate. Then, the wave propeller was turned on, and 
observations concerning the boundary layer were made. This continued for 
different AOA settings and rotational directions of the wave propeller. 

D. RESULTS AND DISCUSSION 

Figures 5.3 and 5.4 show the best pictures taken in the series of 
experiments. Figure 5.3 is a picture of the flow over the stationary wing with 
the wave propeller turned off. The flow is shown to be slightly separated. 
Figure 5.4 shows the flow over the stationary wing with the wave propeller on. 
Here, the flow is seen to be more attached. Although many pictures were 
taken, very few came out clearly because of the poor lighting conditions and 
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limited time available to set appropriate flow visualization conditions for this 
experiment. 

Visual observations of the wave propeller's influence on the stationary 
wing were more convincing than the pictures. It was clear that the wave 
propeller had a beneficial effect on the stationary wing that energized the 
boundary layer at high angle-of-attack and delayed separation. It was observed 
that the rotation direction, AOA, and relative position of the wave propeller all 
controlled the degree of influence on the stationary wing. It was noted that a 
clockwise direction, and increase in AOA had a stronger positive influence on 
the stationary wing. It was unclear as to the best relative position of the wave 
propeller. This was complicated by the limited area available in the test section 
to move the wave propeller. It was further noted that the beneficial influence 
was stronger at lower tunnel speeds. This was a result of the low rotational 
velocities achieved by the wave propeller. Apparently, the ratio of oscillation 
frequency to forward speed is critical for upstream influence. The higher the 
rotational speed of the wave propeller, the faster the freestream velocity can 
be and still achieve a beneficial upstream influence. 
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Figure 5.1 Top view of stationary wing and wave propeller. 




Figure 5.2 Side view of stationary wing and wave propeller. 
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Figure 5.4 Wave propeller is on. 
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VI. CONCLUSION AND RECOMMENDATIONS 



A. SINGLE AIRFOIL ANALYSIS 

The non-linear theory presented here for harmonic motion, and the phase 
lag relationships that exist between the airfoil motion and build-up of 
aerodynamic forces has been extensively verified with existing linear theory. 
Furthermore, a non-linear theory has been presented that predicts the 
propulsive forces and efficiencies associated with pitching, and plunging 
airfoils. Finally, a theory has been presented that suggests a definite link 
between quasi-steady ramp motion (dynamic lift), and a delayed steady state 
event. 

More extensive studies of the quasi-steady phenomenon are recommended 
to thoroughly understand the relationships that exist between steady and 
unsteady motion. Another area of research is the study of the effect of airfoil 
geometry on aerodynamic flutter. A computational method can be developed 
using the presented codes and Theodorsen theory to solve the two-dimensional 
flutter determinant. 

B. TWO AIRFOIL ANALYSIS 

The modified version of USPOTF2A has been verified against existing 
theoretical studies. Furthermore, a limited theory has been presented that 
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predicts the thrust and efficiency associated with two airfoils in close 
proximity. Finally, it was shown in the upwind influence study that an 
oscillating airfoil will have a beneficial influence on the pressure distribution and 
pressure gradient of a stationary airfoil. 

In Chapter III the current limitations of the USPOTF2A code were noted. 
The limitation exists when the discrete vortices from one airfoil's wake interact 
with the other airfoil. The code should be improved to eliminate this difficulty. 
This would greatly enhance the capability of this code when working with large 
amplitude motions or when the airfoils are in tandem. 

C. FLOW VISUALIZATION EXPERIMENTS 

The flow visualization experiment successfully showed the development 
of thrust produced by the wave propeller. Furthermore, the vortex street was 
analyzed to determine which conditions were more favorable for thrust 
production and which conditions induced stall. The enhanced lift flow 
visualization study was not a complete success. The wave propeller did not 
operate smoothly in the vertical position, and the upwind influence was only 
minor. Visual observation indicated a beneficial effect, but further photography 
is needed to document the event. Unfortunately, there was not sufficient time 
for experiments to produce adequate pictures of the phenomenon. 

It is recommended that further experiments be conducted in the low speed 
smoke tunnel. The wave propeller should be modified for a smaller rotational 
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diameter (smaller plunge amplitude). Smaller amplitudes will prevent vortex 
interactions with the tunnel walls. Finally, it is recommended that the wave 
propeller be rebuilt so binding in the bearings will not inhibit the motion in the 
vertical position. Further experiments with smoke visualization and pressure 
measurements are required to better understand the influence of the wave 
propeller on the stationary wing. Certainly, the ability to move the wave 
propeller to various positions around the stationary wing must be included in 
follow-on experiments. 



98 



APPENDIX A 



cccc 

C THIS PROGRAM TAKES THE INPUT FILE (FILE CODE 14) CREATED BY 
C U2DIIF AND CONVERTS THE DATA TO A FREQUENCY, AMPLITUDE, AND 
C PHASE SHIFT. IT ALSO DETERMINES THE PROPULSIVE FORCES AND 
C EFFICIENCY. 

CCC 

PROGRAM PHASESHIFT 

DIMENSION PHASE ( 3 ) ,CL(400) ,CM(400) ,CD(400) , 

+ ALPHA ( 400 ) , TIME (400) ,T(400) ,AMP(3) ,CKT(400) , 

+ FN (400) , R ( 400 ) ,DAT(400) ,FNT(3,400) ,HY(400) 

LOGICAL FLAG 

REAL LI , L2 , L3 , L4 , Ml , M2 , M3 , M4 
PI = ACOS (-1.0) 

C READ TYPE OF MOTION (0=PITCH, 1=PLUNGE) 

READ (14,*) MOTION, ALP1 
C READ NUMBER OF DATA POINTS AND FREQ 
READ (14,*) NPTS , W 
C READ DATA (TIME, CL, CM, CD) 

DO 1=1, NPTS 

READ (8,*) TIME ( I ) , CL ( I ) , CM ( I ) , CD ( I ) 

END DO 

DO 1=1, NPTS 

READ (12,*) T (I) , ALPHA (I) 

END DO 
DO 1=1, NPTS 

READ ( 13 , * ) T ( I ) , HY ( I ) 

END DO 

DO 200, J = 1,2 
DO I = 1 , NPTS 

IF (J .EQ. 1) THEN 
DAT (I) = CL ( I ) 

ELSE IF ( J .EQ. 2) THEN 
DAT ( I ) = CM (I) 

END IF 
END DO 

C READ POSITION DATA 

IF (MOTION .EQ. 1) THEN 
DO 1=1, NPTS 

ALPHA ( I ) = HY ( I ) 

END DO 

ZERO = .00001 
ELSE 

ZERO = .01 
END IF 
15 N=2 
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DO 1=2 , NPTS 

IF (ALPHA ( I ) . LE . ZERO .AND. ALPHA ( I ) .GE. -ZERO) GOTO 16 
N=N+1 
END DO 
16 M=1 

DO I=N+1 ( NPTS 

IF (ALPHA ( I ) . LE . ZERO . AND . ALPHA ( I ) .GE. -ZERO) GO TO 10 
M = M+l 
END DO 

10 CALL AMPLITUDE (DAT, AMP, NPTS, J) 

PRINT*, ' AMPL = ' , AMP ( J) 

DETERMINE PHASE SHIFT 

25 PHI = 0 

ERR = 10000 
TOL = .01 
CN = 1.0 
ITTER = 500 
COUNT = 0 
NUM = 1 

PRINT*, 'N,M,J = ' , N, M, J 

BEGIN ITTERATION TO CONVERGENCE 

30 IF (ERR .LT. TOL) GO TO 35 

SUM = 0 
DO I = N, N+M 

FN(I) = AMP ( J) *SIN (W*T (I) + PHI) 

R ( I ) = ABS ( FN ( I ) - DAT (I) ) 

SUM = SUM + R ( I ) 

END DO 

IF ((SUM - ERR) .GT. 0) THEN 
CN = CN/2.0 
NUM = NUM +1 

IF (NUM/ 2 . 0 .EQ. INT (NUM/2 . 0 ) ) THEN 
FLAG = .FALSE. 

ELSE 

FLAG = .TRUE. 

END IF 
END IF 

IF (FLAG) THEN 

PHI = PHI + CN*PI/180 . 0 
COUNT = COUNT + 1 
ELSE 

PHI = PHI - CN*PI/180 . 0 
COUNT = COUNT + 1 
END IF 

IF (COUNT .GT. ITTER) GOTO 500 
ERR = SUM 
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GOTO 30 

35 PHASE (J) = PHI*180 . 0/PI 
PRINT* , ' ERR = ' , ERR 
GO TO 100 

500 PRINT*, ' ITTERATION GREATER THEN 500. . .INCREASING TOL' 
TOL = 10*TOL 
COUNT = 0 
GO TO 30 

100 IF ( J .EQ. 1) WRITE (6,20) PHASE (J) +180 . 0 

IF ( J .EQ. 2) WRITE (6,21) PHASE(J) 

20 FORMAT ( IX , ' CL PHASE =',F12.6) 

21 FORMAT (IX,' CM PHASE =',F12.6) 

DO I = 4 , NPTS 

FNT ( J , I ) = AMP ( J) *SIN (W*T ( I) + PHASE (J) *PI/l80 . 0 ) 
END DO 
200 CONTINUE 

DO I = 4, NPTS 

WRITE ( 15 , * ) T ( I ) , FNT (1,1) , FNT (2,1) , FNT(3,I) 

END DO 
C 

C DETERMINE THE PROPULSIVE EFFICIENCY 
C 

PHASE (1) = PHASE (1) *PI/180. 0 
PHASE (2 ) = PHASE ( 2 ) *PI/ 180.0 
CDTOT = 0 
K = 0 

DO I =2 , N+M 

CDTOT = CDTOT + CD (I) - CD(1) 

K = K+l 
END DO 

DBAR = CDTOT/ K 
DBAR = DBAR 

PRINT* ,' AVERAGE DRAG (THRUST) , TOTAL DRAG (THRUST) 

+ =', DBAR, CDTOT 

IF (MOTION .EQ. 1) THEN 

WBAR = - .5*W*SIN (PHASE (1) ) *AMP(1) 

ETAS = 2*DBAR/WBAR 
ELSE 

WBAR = . 5 *W*SIN ( PHASE (2 ) ) *AMP ( 2 ) 

ETAS = DBAR /WBAR 
END IF 

PRINT* , ' ETAS , WBAR = ' , ETAS , WBAR 
C 

SUBROUTINE AMPLITUDE (DAT, AMP, NPTS , J) 

DIMENSION DAT (200) , AMP ( 3 ) , AMP1 (10) , AMP2 (10) 

N2 = 0 
M2 = 0 
DO 1= 1,10 

AMPl(I) = 0 
AMP2 (I) =0 
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10 



20 



30 - 



END DO 

DO I = 2 , NPTS - 1 

IF (DAT ( 1+1 ) .LT. 0 .AND. DAT (I) .LT. 0 ) THEN 
IF (ABS (DAT ( 1+1 ) ) .GT.ABS (DAT (I) ) ) THEN 
IF( (N2 + D/2.0 .EQ. INT( (N2 + D/2.0) ) N2 = N2 + 1 
TMP = ABS ( DAT ( I + 1 ) ) 

ELSE 

GOTO 10 
END IF 

IF (TMP .GT. AMP1 (N2 ) ) THEN 
AMP1(N2) = TMP 
TMP =0 
END IF 



ELSE 

IF( (N2 + D/2.0 .NE. INT( (N2 + l)/2.0) ) N2 = 
IF (ABS (DAT ( 1+1 ) ) .GT.ABS (DAT (I) ) ) THEN 
TMP = ABS ( DAT ( I + 1 ) ) 

ELSE 

GOTO 10 
END IF 

IF (TMP .GT. AMP2 (N2 ) ) THEN 
AMP2 ( N2 ) = TMP 

TMP =0 
END IF 
END IF 
END DO 

IF (AMP1 (2) .GT. AMP2 ( 2 ) ) THEN 
IF ( AMP1 ( 2 ) .LT. AMP2 ( 3 ) ) THEN 
COMP = AMP1 (2 ) 

ELSE 



COMP = AMP 2 ( 3 ) 

END IF 
ELSE 

IF (AMP2 (2 ) .LT. AMP1(3)) THEN 
COMP = AMP 2 (2) 

ELSE 



COMP = AMP1 ( 3 ) 

END IF 
END IF 

PRINT* , ' COMP = ' , COMP 
DO I = 2 , N2 

IF(ABS(AMP1(I) -COMP) .GT. .l*COMP) GO TO 20 
M2 = M2 + 1 

AMP ( J) = AMP ( J) +AMP1 (I) 

GO TO 30 

IF (ABS (AMP2 (I) -COMP) .GT. .l*COMP) GO TO 30 
M2 = M2 +1 

AMP ( J) = AMP ( J) + AMP2 (I) 

PRINT*, 'AMP (I) , AMP (I) = AMPl(I), AMP2(I),M2 
END DO 



N2 + 1 
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IF (J .LT. 


3 ) THEN 


AMP (J) 


= AMP ( J) / (M2) 


ELSE 




AMP ( J) 


= AMP (J) / ( 2 *M2 ) 


END IF 




RETURN 




END 
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APPENDIX B 



ccccc 

C THIS PROGRAM TAKES THE INPUT FILE (FILE CODE 14) CREATED BY 
C USPOTF2A AND CONVERTS THE DATA TO A FREQUENCY, AMPLITUDE, 
C AND PHASE SHIFT. IT ALSO COMPUTES THE PROPULSIVE FORCES AND 
C EFFICIENCIES ASSOCIATED WITH THAT DATA SET. 

CCCC 

PROGRAM PHASESHIFT2 

DIMENSION PHASE (2,2) , CL (2 , 450 ) , CM (2 , 450 ) , CD ( 2 , 450 ) , 

+ K(2) , ALPHA (2, 450) , TIME (450) ,T(450) , AMP (2, 2) ,ALP(2) , 

+ CDTOT ( 2 ) , FN (400) ,R(450) ,W(2) ,DAT(2,450) , FNT ( 2 , 450 , 2 ) , 
+ HY (2,450) ,T1 (450) , CPI (450) , CP2 (450) , CP3 (450) , CP4 (450) 
REAL LI , L2 , L3 , L4 , Ml , M2 , M3 , M4 
LOGICAL FLAG, FIRST 
PI = ACOS (-1.0) 

C READ TYPE OF MOTION (0=PITCH, 1=PLUNGE) 

READ (14,*) MOTION, ALP1 
C READ NUMBER OF DATA POINTS AND FREQ 
READ ( 14 , * ) NPTS , W(1),W(2) 

PRINT*, MOTION, ALP1 
PRINT* , NPTS , W ( 1 ) , W ( 2 ) 

C READ DATA (TIME, CL, CM, CD) 

READ ( 7 , * ) 

READ ( 7 , * ) 

READ ( 7 , * ) 

DO 1=1, NPTS 

READ (7, *) TIME (I) , CL (1,1) , CL (2,1) , CM ( 1 , I ) , CM (2 , I) , 

+ CD (1,1) , CD (2,1) 

END DO 
DO 1=1, NPTS 

READ ( 12 , * ) T(I) , ALPHA (1,1) , ALPHA (2, I) 

END DO 
DO 1=1, NPTS 

READ ( 13 , * ) T ( I ) , HY (1,1) , HY (2,1) 

END DO 

IF (MOTION .EQ. 1) THEN 
DO L = 1,2 

DO 1=1, NPTS 

ALPHA(L, I) = HY ( L , I ) 

END DO 
END DO 

ZERO = .00001 
ELSE 

ZERO = .01 
END IF 

IF (W(l) .EQ. 0) THEN 



104 



FIRST = .FALSE. 

W(l) = W(2) 

DO I = 1 , NPTS 

ALPHA (1,1) = ALPHA (2,1) 

END DO 

ELSE IF (W(2) .EQ. 0) THEN 
FIRST = . TRUE . 

W(2) = W(l) 

DO I = 1, NPTS 

ALPHA (2,1) = ALPHA (1,1) 

END DO 
END IF 

DO 200, J = 1,2 
DO I = 1 , NPTS 

IF (J .EQ. 1) THEN 

DAT (1,1) = CL (1,1) 

DAT (2, I) = CL (2, I) 

ELSE IF ( J .EQ. 2) THEN 
DAT (1,1) = CM (1,1) 

DAT (2,1) = CM (2,1) 

END IF 
END DO 

C READ POSITION DATA 
DO 100 L = 1,2 

15 N=2 

DO 1=2, NPTS 

IF (ALPHA(L, I) .LE. ZERO .AND. ALPHA (L, I) .GE. 

+ -ZERO) GO TO 16 

N=N+1 
END DO 

16 M=1 

DO I=N+1 , NPTS 

IF (ALPHA (L, I) .LE. ZERO .AND. ALPHA (L, I) .GE. -ZERO) 
+ GO TO 10 
M = M+l 
END DO 

10 PRINT*, W(L) 

CALL AMPLITUDE (DAT, AMP, NPTS, L,J,N) 

PRINT*, ' AMPL = ' , AMP ( J , L ) 

C 

C DETERMINE PHASE SHIFT 
C 

25 PHI = 0 

ERR = 10000 
TOL = .01 
CN = 1.0 
ITTER = 500 
COUNT = 0 
NUM = 1 
FLAG = .TRUE. 
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PRINT*, 'N,M,L,J = ' , N,M, L, J 
C BEGIN ITTERATION TO CONVERGENCE 
30 IF (ERR .LT. TOL) GO TO 35 

SUM = 0 
DO I = N, N+M 

FN ( I ) = AMP ( J, L) *SIN (W (L) *T ( I ) + PHI) 

R ( I ) = ABS(FN(I) - DAT ( L , I ) ) 

SUM = SUM + R ( I ) 

END DO 

IF ( (SUM - ERR) .GT. 0.0) THEN 
CN = CN/2.0 
NUM = NUM + 1 

IF (NUM/ 2 . 0 .EQ. INT (NUM/ 2.0)) THEN 
FLAG = .FALSE. 

ELSE 

FLAG = .TRUE. 

END IF 
END IF 

IF (FLAG) THEN 

PHI = PHI + CN*PI/180 . 0 
COUNT = COUNT + 1 
ELSE 

PHI = PHI - CN*PI/180 . 0 
COUNT = COUNT + 1 
END IF 

IF (COUNT .GT. ITTER) GOTO 500 
ERR = SUM 
GOTO 30 

35 PHASE ( J, L) = PHI*180/PI 

PRINT*, ' PHASE (J,L) = ', PHASE (J,L) 

PRINT* , ' ERROR = ' , ERR 

GO TO 100 

500 PRINT* , ' ITTERATION GREATER THEN 500 .. . INCREASING TOL' 
TOL = 10*TOL 
CN = 10*CN 
COUNT = 0 
GO TO 30 
100 CONTINUE 

IF ( J .EQ. 1) WRITE (6,20) PHASE (J, 1) +180 . 0 ,PHASE(J,2) 





+ 


+ 180.0 

IF ( J .EQ. 2) 


WRITE (6,21) 


PHASE ( J, 1) , 


PHASE ( J, 2 ) 


20 


+ 


FORMAT ( IX , ' 
=' , F10 . 6 ) 


CL 


PHASE (1) = 


' , F10 . 6 , 5X, ' 


CL PHASE (2) 


21 


+ 


FORMAT (IX, ' 
=' , F10 . 6 ) 


CM 


PHASE (1) = 


' , F10 . 6 , 5X, ' 


CM PHASE (2) 



DO I = 4 , NPTS 

FNT ( J, 1,1) =AMP ( J, 1) *SIN(W(1) *T(I) + PHASE (J, 1) *PI/180.0) 
FNT ( J , I , 2 ) = AMP ( J , 2 ) *SIN(W(2) *T(I) + PHASE ( J, 2 ) *PI/180 . 0) 
END DO 

PHASE ( J, 1) = PHASE (J, 1) *PI/180.0 
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PHASE ( J, 2 ) = PHASE (J,2) *PI/180.0 
200 CONTINUE 

DETERMINE THE PROPULSIVE EFFICIENCY 
INPUT DTS 

PRINT* , ' INPUT DTS ' 

READ* , DTS 
K1 = 0 
DO L = 1,2 
CDTOT(L) = 0 
K (L) = 0 

DO I =DTS , 2*DTS 

CDTOT(L) = CDTOT(L) + CD(L,I) 

K (L) = K (L) + 1 
END DO 
END DO 

DBAR1 = CDTOT ( 1 ) / K ( 1 ) 

DBAR2 = CDTOT ( 2 ) / K ( 2 ) 

PRINT* AVERAGE DRAG , TOTAL DRAG (AF 1) 
+ ' ,DBAR1, CDTOT (1) , K(l) 

PRINT* AVERAGE DRAG, TOTAL DRAG (AF 2) 
+ ' , DBAR2 , CDTOT ( 2 ) ,K(2) 

IF (MOTION .EQ. 1) THEN 
IF (FIRST) THEN 
PHAS = PHASE (1,1) 

AMPL = AMP (1,1) 

ELSE 

PHAS = PHASE (1,2) 

AMPL = AMP (1,2) 

END IF 

WBAR = - . 5 * W ( 1 ) * S IN ( PHAS ) * AMPL 
ELSE 

IF (FIRST) THEN 
PHAS = PHASE (2,1) 

AMPL = AMP (2,1) 

ELSE 

PHAS = PHASE (2, 2) 

AMPL = AMP (2,2) 

END IF 

ALP1 = ALP1*PI/ 180.0 
WBAR = - . 5*W ( 1) *SIN (PHAS) *AMPL 
END IF 

ETAS = 2* (DBAR1+DBAR2) /WBAR 
PRINT* , ' ETAS , WBAR = ' , ETAS , WBAR 
DO I=DTS , 2*DTS 

CL1TOT = CL1TOT + CL (1,1) 

CL2TOT = CL2TOT + CL (2, I) 

K1 = K1 + 1 
END DO 
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CL1AVG = CL1T0T/K1 
CL2AVG = CL2TOT/K1 
DELCL1 = CL1AVG - CL (1,1) 

DELCL2 = CL2AVG - CL (2,1) 

PRINT*, 'CL 1AVG, DELCL1 = ' , CL1AVG , DELCL1 
PRINT*, 'CL2AVG, DELCL2 = ' , CL2AVG, DELCL2 
C 

DO I = 4 , NPTS 

WRITE (15, *)T(I) ,FNT( 1,1,1) ,FNT(1,I,2) ,FNT(2,I,1) , 
+ FNT (2,1,2) 

END DO 
C 

END 
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SUBROUTINE AMPLITUDE (DAT, AMP, NPTS , L, J, N) 

DIMENSION DAT (2,200) ,AMP(2,2) , AMP1 (10) , AMP2 (10) 

N2 = 0 
M2 = 0 
DO 1= 1,10 

AMPl(I) = 0 
AMP 2 (I) =0 
END DO 

DO I = 2, NPTS -1 

I F ( DAT ( L , I + 1 ) .LT. 0 .AND. DAT ( L , I ) .LT. 0 ) THEN 
IF (ABS (DAT (L, 1+1) ) .GT. ABS (DAT (L, I) ) ) THEN 
IF( (N2+l)/2.0 .EQ. INT( (N2+l)/2.0) ) N2 = N2+1 
TMP = ABS ( DAT ( L , I + 1 ) ) 

ELSE 

GOTO 10 
END IF 

IF (TMP .GT. AMP1 (N2 ) ) THEN 
AMP1(N2) = TMP 
TMP =0 
END IF 



ELSE 

IF( (N2 + l)/2.0 .NE. INT( (N2 + D/2.0) ) N2 =N2 + 1 
IF (ABS (DAT(L, 1+1) ) . GT . ABS (DAT ( L, I ) ) ) THEN 

TMP = ABS ( DAT ( L , I + 1 ) ) 

ELSE 

GOTO 10 
END IF 

IF (TMP .GT. AMP2 (N2 ) ) THEN 
AMP2 (N2 ) = TMP 
TMP =0 
END IF 
END IF 
END DO 

IF (AMP1 (2) .GT. AMP2 ( 2 ) ) THEN 
IF ( AMP1 ( 2 ) .LT. AMP2 ( 3 ) ) THEN 
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COMP = AMP1 (2) 

ELSE 

COMP = AMP2 ( 3 ) 

END IF 
ELSE 

IF ( AMP2 ( 2 ) .LT. AMP1(3)) THEN 
COMP = AMP2 ( 2 ) 

ELSE 

COMP = AMP1 (3) 

END IF 
END IF 

PRINT* , ' COMP = ' , COMP 
DO I = 2 , N2 

IF (AMP1 ( I ) .GT. COMP) GO TO 19 

IF(ABS (AMP1 (I) -COMP) .GT. .12*COMP) GO TO 20 



19 


M2 = M2 + 1 






AMP ( J, L) = AMP ( J, L) 
GO TO 30 


+AMP1 (I) 


20 


IF (AMP2 ( I ) .GT. COMP) GO 


TO 21 




IF (ABS ( AMP2 (I) -COMP) 


.GT. . 12 *COMP) GO TO 30 


21 


M2 = M2 +1 






AMP ( J, L) = AMP ( J , L ) 


+ AMP 2 (I) 


30 


PRINT* , ' AMP ( I ) , AMP ( I ) 
END DO 

AMP ( J , L ) = AMP ( J , L ) / ( M2 ) 

RETURN 

END 


= AMPl(I), AMP2 ( I ) , M2 
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